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Directed and elliptic flow measurements for charged pions and protons are reported as a function of 
transverse momentum, rapidity, and centrality for 40 and 158^4 GeV Pb + Pb collisions as recorded 
by the NA49 detector. Both the standard method of correlating particles with an event plane, and 
the cumulant method of studying multiparticle correlations are used. In the standard method the 
directed flow is corrected for conservation of momentum. In the cumulant method elliptic flow is 
reconstructed from genuine 4, 6, and 8-particle correlations, showing the first unequivocal evidence 
for collective motion in A+A collisions at SPS energies. 
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I. INTRODUCTION 

In non-central collisions, collective flow leads to characteristic azimuthal correlations between particle momenta 
and the reaction plane. The geometry of a non-central collision between spherical nuclei is uniquely specified by the 
collision axis and the impact parameter vector b. In particular, the latter defines a unique reference direction in 
the transverse plane. Directed flow (vj) and elliptic flow (1/2) cause correlations between the momenta of outgoing 
particles with this reference direction [lj. They are defined by 

Vi = (cos(0 - $ RP )) , v 2 = (cos(2(</> - $rf))) , (I) 

where 4> denotes the azimuthal angle of an outgoing particle and Qrp is the azimuthal angle of b. Angular brackets 
denote an average over particles and events. 

Since primary collisions between the incoming nucleons are expected to be insensitive to the direction of impact 
parameter at high energy, azimuthal correlations are believed to result from secondary interactions, or final state 
interactions. As such, they are sensitive probes of "thermalization," which should be achieved if final state interactions 
are strong enough. In the energy regime in which relatively few new particles are created, the flow effects are due to 
the nucleons which participate in the collision. Thus at low energy flow is used to study the properties of compressed 
nuclear matter and more specifically the nuclear equation of state |2( . In nuclear collisions at ultrarelativistic energies 
the number of newly created particles is so large that their behavior will dominate the observable flow effects. 

For the interpretation of experimental results on flow, theoretical tools are needed. There are two types of models to 
describe final state interactions, based either on hydrodynamics or on a microscopic transport (or cascade) approach. 
Hydrodynamics is adequate when the mean free path of particles is much smaller than the system size. Then, 
the interactions between the various particles in the system can be expressed in terms of global thermodynamic 
quantities, i.e., an equation of state. In this essentially macroscopic description, the collective motion results from a 
pressure gradient in the reaction volume, the magnitude of which depends upon the compressibility of the underlying 
equation of state Since partonic and hadronic matter are expected to have different compressibilities, it may 
be possible to deduce from a flow measurement whether it originates from partonic or hadronic matter, or from the 
hadronization process taking place during the transition between the two [2 S S @ ■ Microscopic cascade models 
are more appropriate when the mean free path of the particles is of the same order, or larger, than the size of the 
system, which is often the case in heavy ion collisions. They require a more detailed knowledge of the interactions 
(cross sections, etc.) of the various particles in the medium. Detailed flow analyses may help to falsify or confirm the 
corresponding model assumptions. 

In particular, the study of the energy dependence of flow is considered 0, H, as a promising strategy in the 
search for evidence for the hypothesis that the onset of deconfinement occurs at low SPS energies. More generally, 
the energy scan project [nj at the CERN Super Proton Synchrotron (SPS) was dedicated to the search for the onset 
of deconfinement in heavy ion collisions. In fact anomalies observed in the energy dependence of total kaon and pion 
yields can be understood as due to the creation of a transient state of deconfined matter at energies larger than 
about 40A GeV 

In that context, this paper presents the most detailed analysis so far of directed and elliptic flow of pions and 
protons at various SPS energies. The first publication from NA49 on anisotropic flow [l3| was based on a small set of 
I58A GeV data with a medium impact parameter trigger with only the tracks in the Main TPCs used in the analysis. 
Subsequently a method was found for improving the second harmonic event plane resolution and revised results were 
posted on the web 14]. The present analysis of 1 58 A GeV data is both more detailed and more accurate than the 
previous one: it uses much larger event statistics, a minimum bias trigger, integration over transverse momentum p t 
or rapidity y using cross sections as weights, and improved methods of analysis. Moreover, in this paper the NA49 
results on flow at 40A GeV are published for the first time. Preliminary results from this analysis have been presented 
in Refs. [I1GIG3. 

Two types of methods are used in the flow analysis: the so-called "standard" method 0,0, H(| requires for each 
individual collision an "event plane" which is an estimator of its reaction plane. Outgoing particles are then correlated 
with this event plane. This method, however, neglects other sources of correlations: Bose-Einstein (Fermi-Dirac) 
statistics, global momentum conservation, resonance decays, jets, etc. The effects of these "nonflow" correlations may 
be large at the SPS, as shown in Ref. |2lL l22l|. The standard method has been improved to take into account part 
of these effects. In particular, correlations from momentum conservation are now subtracted following the procedure 
described in Ref. [23|]. Recently, a new method has been proposed which allows to get rid of nonflow correlations 
systematically, independent of their physical origin [241 I2R l2q . This method extracts directed and elliptic flow 
from genuine multiparticle azimuthal correlations, which are obtained through a cumulant expansion of measured 
multiparticle correlations. The results obtained with both methods will be presented and compared. 

The paper is organized as follows: Sec. ITT1 covers the experiment, Sec. IIIII describes the data sets, the selection 
criteria for events and particles, and the acceptance of the detector, in Sec. lIVI the two methods of flow determination 
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are explained, and Sec. contains the results on elliptic and directed flow as function of centrality and beam energy. 
Sec. I VII focuses the discussion on model calculations and Sec. IV 111 summarizes the paper. 



II. EXPERIMENT 

The NA49 experimental set-up [2?} is shown in Fig. ^ It consists of four large-volume Time Projection Chambers 
(TPCs). Two of these, the Vertex TPCs (VTPC-1 and VTPC-2), are placed in the magnetic field of two super- 
conducting dipole magnets (VTX-1 and VTX-2). This allows separation of positively and negatively charged tracks 
and a precise measurement of the particle momenta (p) with a resolution of a(p)/p 2 = (0.3 — 7) • 10~ 4 (GeV/c) -1 . The 
other two TPCs (MTPC-L and MTPC-R), positioned downstream of the magnets, were optimized for high precision 
detection of the ionization energy loss dE/dx (relative resolution of about 4%) which provides a means to measure 
the particle mass. The TPC data yield spectra of identified hadrons above midrapidity. The magnet settings at 158A 
GeV were: £ (VTX-1) ^ 1.5 T and 5(VTX-2) = 1.1 T. In order to optimize the NA49 acceptance at AOA GeV the 
magnetic fields of VTX-1 and VTX-2 were lowered in proportion to the beam momentum. Data were taken for both 
field polarities. 
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FIG. 1: (Color online) The experimental set-up of the NA49 experiment |27|. 

The target (T), a thin lead foil (224 mg/cm 2 , approx. 0.47% of Pb-interaction length), was positioned about 80 cm 
upstream from VTPC-1. Beam particles were identified by means of their charge as seen by a gas Cherenkov counter 
(S2') in front of the target. An identical veto-counter directly behind the target (S3) is used to select minimum bias 
collisions by requiring a reduction of the Cherenkov signal by a factor of about 6. Since the Cherenkov signal is 
proportional to Z 2 , this requirement ensures that the projectile has interacted with a minimal constraint on the type 
of interaction. This limits the triggers on non-target interactions to rare beam-gas collisions, the fraction of which 
proved to be small after cuts, even in peripheral Pb+Pb collisions. The counter gas which is present at atmospheric 
pressure all the way from S2' to the target and further on to S3 was He for the 158A GeV and CO2 for the 40A 
GeV runs. The signal from a small angle calorimeter (VCAL), which measured the energy carried by the projectile 
spectators, was used to make off-line centrality selections. The geometrical acceptance of the VCAL calorimeter was 
adjusted for each energy in order to cover the projectile spectator region by a proper setting of a collimator (COLL) 
[271 |28j • The NA49 coordinate system is defined as right handed with the positive z-axis along the beam direction 
and the x-axis in the horizontal and the y-axis in the vertical plane. 
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III. DATA 

The data on Pb+Pb collisions at 40^4 GeV were collected within the energy scan program at the CERN SPS [To| . 
As part of this program Pb+Pb collisions at 20, 30, 40 and 80^4 GeV were recorded by the NA49 detector during 
heavy ion runs in 1999, 2000 and 2002. Due to limited beam time, minimum bias data required for a flow analysis 
were not taken at 80 A GeV, and the 20 and 30 A GeV data from 2002 have not yet been analyzed. The corresponding 
data at the top SPS energy (158^4 GeV) were taken from runs in 1996 and 2000. 

A. Data sets 

The data sets used in this analysis were recorded at 158A GeV and 40A GeV with a minimum bias trigger allowing 
for a study of centrality and energy dependence. Since central collisions have a small weight in such a selection of 
events, their number was augmented by data from central trigger runs at 158A GeV. The final results for the 40A 
GeV beam minimum bias data were obtained from 350 k minimum bias events for the standard method and 310 
k events for the cumulant method. For the 158A GeV results, the minimum bias events used were 410 k for the 
standard method and 280 k for the cumulant method. The 12.5% most central events added in were 130 k for the 
standard method and 670 for the cumulant method. In addition, for the integrated cumulant results, 280 k events 
from another run triggered on 20% most central collisions was added. These numbers refer to events which fulfill the 
selection criteria. After verifying that the analysis of events recorded with opposite field polarities give compatible 
results the corresponding data sets were combined and processed together. Full coverage of the forward hemisphere 
for pions and protons is achieved by using the tracks combined from both the Vertex and Main TPCs. 

B. Selections and particle identification 

The sample of events provided by the hardware trigger is contaminated by non-target interactions which are removed 
by a simultaneous cut on the minimum number of tracks connected to the reconstructed primary vertex (10) and on 
the deviation from its nominal position in space (0.5 cm in all dimensions). Quality criteria ensured that only reliably 
reconstructed tracks were processed. This acceptance was defined by selecting tracks in each TPC if the number of 
potential points in that TPC based on the geometry of the track was at least 20 in the vertex TPCs and 30 in the 
main TPC. In order to avoid split tracks the number of fit points for the whole track had to be greater than 0.55 
times the number of potential points for that track. The x 2 P er degree of freedom of the fit had to be less than 10. 
Tracks with transverse momentum (p t ) up to 2 GeV/c were considered. The fraction of tracks of particles from weak 
decays or other secondary vertices was reduced by cutting on the track distance from the reconstructed event vertex 
in the target plane (±3 cm in the bending and ±0.5 cm in the non-bending direction). 

The binning of the event samples in centrality was done on the basis of the energy measurement in the forward 
calorimeter (VCAL). Its distribution was divided into 6 bins with varying widths. Each bin has a mean energy (E ) 
and corresponds in a Glauber-like picture to an impact parameter range (b) with an appropriate mean, a mean number 
of wounded nucleons {Nwn}, a mean number of participants (N part ), and a cross section fraction ct/ctt with cry being 
the total hadronic inelastic cross section of Pb+Pb collisions which has been estimated to be 7.1 b at both energies. 
Details of the binning are given in Table [I] In the graphs "central" refers to bins 1 plus 2, "mid-central" to bins 3 
plus 4, and "peripheral" to bins 5 plus 6. When we integrate over the first five centrality bins to present "minimum 
bias" results, we believe we have integrated out to impact parameters of about 10 fm corresponding to cr/cry=0.435. 

In the standard method of flow analysis the determination of the event plane is required (see below) . The uncertainty 
of its azimuthal angle in the laboratory coordinate system depends not only on the total number of particles used but 
also on the size and sign of the flow signal of these particles, which are in general different for different types of flow 
and different phase space regions. To optimize the resolution of the event plane orientation the following selection of 
tracks used for the event plane determination was made using y in the center of mass and pt- For the first harmonic: 
< p t < 1 GeV/c (centrality bins 3-6), < p t < 0.3 GeV/c (centrality bin 1), < p t < 0.6 GeV/c (centrality bin 
2), 1.1 < y < 3.1 for 158A GeV data, and < p t < 1 GeV/c (centrality bins 1-6), 0.8 < y < 2.8 for 40A GeV 
data. For the second harmonic: < p t < 1 GeV/c, —0.5 < y < 2.1 for 158A GeV data, and < p t < 1 GeV/c, 
—0.4 < y < 1.8 for 40A GeV data. In comparison to these selections for good event plane resolution, the differential 
data to be presented go to higher p t but lower maximum y values. 

Particle identification is based on energy loss measurements (dE/dx) in the time projection chambers. An enriched 
sample of pions is obtained by removing those particles which are obviously not pions by appropriate cuts in the lab 
momentum - dE/dx plane. The remaining contamination amounts to less than 5% for negatively charged pions. For 
positively charged pions it is less than 20% between 2 and 20 GeV/c momentum in the laboratory. Outside this range 
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the contamination increases up to 35% for lower momenta. For higher momentum the contribution to the measured 
flow is small due to the vanishing cross section of pions in this region. Although the fraction of misidentified particles 
is substantial, the effect on the results will be small, since v-i for pions and protons is comparable and depends on 
rapidity and pt in a similar way. The kaons are expected to follow the same trend. To examine the influence of 
the contamination, v± and V2 of all charged particles were analyzed and compared to results for pions. The small 
differences are included in the quoted systematic errors. The proton identification is restricted to laboratory momenta 
above 3 GeV/c and thus to the region of the relativistic rise of the specific energy loss. Tight upper limits of dE/dx 
remove almost quantitatively all lighter particles. The remaining contamination amounts to less than 5% kaons and 
pions for 158A GeV data and less than 8% for 4QA GeV data. 

TABLE I: Listed for the two beam energies and six centralities are: E° / E beam , the forward calorimeter energy divided by the 
beam energy; (E°/E beam ), the mean value; a/ar, the fraction of the total cross section in that bin; the integral of a/ar; b, the 
estimated range of impact parameters; (6), the estimated mean impact parameter; (Nwn), the estimated number of wounded 
nucleons, and (N par t), the estimated mean number of participants. The 6 bins in centrality have different values of E°/E beam 
at 40 and 158A GeV since the acceptance of the forward calorimeter depends on beam energy. 
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C. Acceptance 

The NA49 detector was designed for large acceptance in the forward hemisphere of the center-of-mass frame. The 
resulting acceptance is illustrated by the density distributions of protons and pions as a function of rapidity and 
transverse momentum as seen in Ref. |ll|. The scaling of the magnetic field strength with beam energy ensures 
similar distributions at both energies. The NA49 detector employs two dipole magnets with main field components 
perpendicular to the beam axis. This breaking of rotational symmetry together with the rectangular TPC shapes 
introduces azimuthal anisotropies which are more pronounced at the lower beam energy. The Lorentz boost focuses 
the tracks of all particles forward of midrapidity into cones of approximately 5 and 10 degrees at 158 and 40^4 GeV, 
respectively. Acceptance losses occur for particles at large angles with respect to the bending plane. A typical inclusive 
azimuthal angle (0) distribution is shown in Fig. [21 top for pions, which in an ideal detector would be flat. 

In the standard method, event plane determinations from such a distribution would obviously lead to acceptance 
biased event planes. As long as each bin in the acceptance is populated significantly, the bias can be removed as 
described below by recentering the particles in the plane perpendicular to the beam for each rapidity and pt bin, in 
such a way that the event plane distribution becomes flat. The result of this procedure is exemplified for the first 
harmonic in Fig. |3 bottom. 

In the cumulant method, the main contributions of the detector inefficiency, that is, spurious correlations which have 
nothing to do with physical (flow or nonflow) correlations, are automatically removed. The only required acceptance 
corrections amount to a global multiplicative factor which depends solely on the specific detector under study, and 
can thus be calculated separately from the flow analysis [2j|. Further details will be given in Sec. lIVBl 
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FIG. 2: (Color online) On top is a typical azimuthal distribution of particles (pions, third centrality bin) at 40A GeV. On the 
bottom is the azimuthal distribution of the first harmonic event planes after correction for the laboratory azimuthal anisotropies. 



IV. METHODS 



In this section, we recall the principles of the standard (Sec. IIV All and cumulant (Sec. IIV Bp methods of flow 
analysis. 



A. Standard method 

The standard method |l9l l20j| correlates the azimuthal angles of particles, 4>, with an estimated event plane to 
obtain the observed coefficients in a Fourier expansion in the plane transverse to the beam. The observed coefficients 
are then divided by the resolution of the event plane obtained from the correlation of the estimated event planes 
of two random subevents. The estimated event plane angle, <!>„, are obtained from the azimuthal angles of the Q n 
vectors whose x and y components are defined by: 

Q n cos n(f>„ = Wi (cos n<f>i — (cos ncj)) ) , 

i 

QnS'mn^n = W i ( Sm n 4>i ~ ( Sm n §) ) i ( 2 ) 

i 

where n is the harmonic order and the sum is taken over the M particles in the event. In this work the weights, 
Wi, have been taken to be pt for the second harmonic and y in the center of mass for the first harmonic. To make 
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the event plane isotropic in the laboratory in order to avoid acceptance correlations, we have used the recentering 
method |20j . The mean (sin ncj)) and (cos n<j)) values in the above equation were calculated as a function of p t and y 
for all particles in all events in a first pass through the data, and then used in a second pass to recenter the Q n vector 
to be isotropic as shown in Fig. [21 bottom. The mean sin and cos values were stored in a matrix of 20 pt values and 
50 y values for each harmonic. Particles were only used for the event plane determination if the absolute value of the 
mean sin and cos values for that bin were less than 0.2. Then the flow values are calculated by 

„. _ (cos(n(& - $„))) ^ 



(cos(n($„ - $i? P ))) ' 



If particle i was used also for the event plane determination, its contribution to Q n is subtracted before calculating 
$„, so as to avoid autocorrelations. The denominator is called the resolution and corrects for the difference between 
the estimated even t plane and the rea l reaction plane, <f>np. It is obtained from the resolution of the subevent event 
planes, which is \J (cos(n($ a — The resolution of the full event plane, for small resolution, is approximately 

y/2 larger, but the actual equation in Ref. [2£j was used for this calculation. It was found to be more accurate to 
calculate relative to the second harmonic event plane, $2, although the sign of vi was determined to be positive by 
correlation with the first harmonic event plane $1. The sign of v\ was set so that protons at high rapidity have positive 
V\ as described below. The software used in this analysis was derived from that used for the STAR experiment p9|. 

Equation @ is the most general form to determine v n . In the case of the NA49 experiment the main losses are 
concentrated around 90 deg and 270 deg in the up and down directions. (See Fig. top.) In order to limit the 
analysis to the regions of more uniform acceptance a cut on the particle azimuthal angle was applied: Particles with 
cos(20) < are cut out. This however requires large acceptance corrections if Eq.[3]is used for V2 determination. The 
correlation term may be modified in Eq. in order to select azimuthal regions with small distortions. Its numerator 
may be rewritten in the form of a sum of products instead of a difference of angles: 

(cos(n0i) • cos(n$„)) + (sin(n<&) • sin(n$„)) . (4) 

Since the (ideal) inclusive azimuthal distributions in A +A collisions are flat by definition, both terms must be of equal 
magnitude and v n can be calculated by 2 times either the first or second term. The cut applied on the azimuthal 
angle (cos(2</>) < 0) leads to large corrections to the cosine term in Eq. J3J. Therefore, for the second harmonic at 
40 A GeV the numerator in the vi calculation was done using the expression 

2 (sin(n0i) • sin(n$„)) . (5) 

The same line of argument applies to the denominator of Eq. J2J where again only the sin terms were used and the 
subevent plane resolution was increased by a factor of y/2: 

_ 2 • (sin(n0i) • sin(n$„)) 



V2 • (sin(n$„) • sm(n$ RP )) 

Determination of V2 by Eq. requires the acceptance correction only for losses in the selected ang ular ranges. 

A momentum conservation correction for the first harmonic was made as described in Ref. [23| . The correction is 
made to the observed differential flow values before they are divided by the event plane resolution. Figure shows 
that this correction, without any adjustable parameter, makes the directed flow curve cross zero at midrapidity. The 
figure also shows that because this correction is proportional to pt there is a large effect at high p t . There is no effect 
on the elliptic flow because it is calculated relative to the second harmonic event plane. Table[n]shows the parameters 
used in making this correction. 



B. Cumulant method 



In this section we first recall the motivations for developing alternative methods of flow analysis. We then explain 
the principle of the method. Unlike the standard method, the cumulant method yields in principle several independent 
estimates of directed and elliptic flow, which will be defined below. Finally, we describe the practical implementation 
of the method. 

At the core of the standard method outlined in Sec. IIV Al lies a study of two-particle correlations: one correlates 
either two subevents, to derive the event plane resolution, or one particle with (a second particle belonging to) the Q n 
vector. The basic assumption is that the correlation between two arbitrary particles is mainly due to the correlation 
of each single particle with the reaction plane, that is, due to flow. However, there exist other sources of two-body 
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FIG. 3: (Color online) Charged pion directed flow as a function of rapidity (top) and pt (bottom) for minimum bias 158A 
GeV Pb + Pb. Shown are vi before (squares) and after (circles) correction for momentum conservation. The open points in 
the top graph have been reflected about midrapidity. The solid lines are polynomial fits. 

correlations, that do not depend on the reaction plane; for instance, physical correlations arising from quantum (HBT) 
effects, global momentum conservation, resonance decays, or jets. At SPS energies, it turns out that these "nonflow" 
two- particle correlations are a priori of the same magnitude as the correlations due to flow [2lL |22| , at least in some 
phase space regions. While some of these correlations can be taken into account in the standard method with a 
minimal modeling of the collisions (see end of Sec. If VAT) , others cannot be estimated as reliably. 

This observation motivated the elaboration of new methods of flow analysis, which are much less biased by nonflow 
correlations than the standard method [24|. The basic idea of the methods is to extract flow from multiparticle 
azimuthal correlations, instead of using the correlation between two particles only. Naturally, the measured fc-body 
correlations also consist of contributions due to flow and to nonflow effects. Nevertheless, by performing a cumulant 
expansion of the measured correlations, it is possible to disentangle the flow contribution from the other, unwanted 
sources of correlations. Thus, at the level of four-particle correlations, one can remove all nonflow two- and three- 
particle correlations, keeping only the correlation due to flow, plus a systematic uncertainty arising from genuine 
nonflow /owr-particle correlations, which is expected to be small. 

The cumulant method not only minimizes the influence of nonflow correlations; it also provides several independent 
estimates of v\ and «2, which will be labeled by the order k at which the cumulant expansion is performed: for instance, 
«2{4} denotes our estimate of «2 using cumulants of 4-particle correlations, etc. Generally speaking, the systematic 
error due to nonflow correlations decreases as the order k increases, at the expense of an increased statistical error. 

To be more specific, we first consider a simplified situation where one wishes to measure the average value of the flow 
v n over the detector acceptance, which is assumed to have perfect azimuthal symmetry. The lowest order estimate of 




p, < 2 GeV/c 
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TABLE II: Listed for the two beam energies and six centralities are: N, the estimated total multiplicity of charged plus neutral 
hadrons over all phase space; (pt), the estimated mean p 2 of all hadrons over all phase space; M, the mean multiplicity of 
particles used for the event plane determination; /, the fraction defined by Eq. (17) in Ref. 23], which controls the correction 
for momentum conservation; Xi, the first harmonic resolution parameter for the full event plane (it follows the convention of 
Ref. {Ifflj and is y/2 smaller than defined in Ref. [2fT|'l: the resolution of the first harmonic full event plane; the percent increase in 
the resolution due to momentum conservation, and the resolution of the second harmonic event plane. The first two centralities 
at 158A GeV had more restrictive p t cuts. The first centrality bin at 40A GeV did not have sufficient statistics to determine 
the second harmonic event plane. 



Centrality 


1 


2 


3 


4 


5 


6 


158A GeV 














N 


2402 


1971 


1471 


1028 


717 


457 


(pi) (GeV 2 /c 2 ) 


0.32 


0.32 


0.31 


0.30 


0.28 


0.27 


M 


119 


181 


154 


110 


78 


46 


f 


0.07 


0.14 


0.17 


0.17 


0.18 


0.17 


XI 


0.25 


0.27 


0.40 


0.44 


0.48 


0.47 


resolution 1st 


0.22 


0.24 


0.34 


0.38 


0.40 


0.40 


res. increase, % 


4 


16 


11 


9 


7 


7 


resolution 2nd 


0.23 


0.28 


0.36 


0.40 


0.37 


0.29 


40A GeV 














N 


1473 


1215 


913 


643 


453 


290 


(p?) (GeV 2 /c 2 ) 


0.39 


0.38 


0.38 


0.36 


0.34 


0.32 


M 


89.3 


76.4 


59.5 


42.2 


29.7 


17.2 


f 


0.11 


0.11 


0.11 


0.12 


0.12 


0.11 


Xi 


0.23 


0.23 


0.30 


0.34 


0.38 


0.40 


resolution 1st 


0.20 


0.20 


0.26 


0.29 


0.33 


0.34 


res. increase, % 


14 


16 


8 


6 


5 


4 


resolution 2nd 




0.13 


0.25 


0.31 


0.25 


0.21 



v n from two-particle correlations, v n {2}, is then defined by 

v n {2} 2 = (e'*"^)) (7) 



where brackets denote an average value over pairs of particles emitted in a collision, and over events. Please note that 
v n {2} is a priori consistent with the value given by the standard method, Eq. @, at least if the cuts in phase space 
are identical in both analyses. 

Higher order estimates are obtained from two complementary multiparticle methods. The first one [25] measures 
the flow harmonics separately, either v% or V2- For instance, the four-particle estimate t> n {4} is defined by 



v n {4} 4 



m(0i+0 2 — 03— <t>i) 



e in(0i-0 3 )\ / e in(rj> 2 -<t>4,) 



e in (*»-*0\, (8) 

where the average runs over quadruplets of particles emitted in the collision, and over events. The right-hand side 
defines the cumulant of the four-particle correlations. This can be generalized to an arbitrary even number of particles, 
which yields higher order estimates v„\p\, i> n {8}, etc. 

The second multiparticle method |26( was used to analyze directed flow v\. It relies on a study of three-particle 
correlations which involve both v± and v^: 



e 



H*+**-*M) ~ (9) 



In the case of NA49, we shall see that the first multiparticle method provides reliable estimates of v-i (the most reliable 
was found to be t>2{2} at 40A GeV and «2{4} at 158A GeV, as will be discussed later). Then, the above equation can 
be used to obtain an estimate of v±, which is denoted by ^i{3} since it involves a 3-particle correlation. As shown 
in Ref. [2(j, and will be seen below in Sec. IV Bl fi{3} offers the best compromise between statistical errors (which 
prevent obtaining i>i{4} with the first method) and systematic errors from nonflow correlations, which plague the 
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lowest order estimate ^i{2}. In particular, among other nonflow correlations, fi{3} is insensitive to the correlation 
due to momentum conservation, so that one need not compute it explicitly as in the case of the standard method. 
As a matter of fact, a straightforward calculation using the i/iree-particle correlation due to momentum conservation, 
given by Eq. (12) in Ref. j^l), shows that the contribution of pt conservation to the average ^e l( -^ 1+ ^ 2 ~ 2 ^ 3 )) is 

(Pt) 2 /(N 2 (Pt )), roughly smaller than l/N 2 . (Please note that this three-particle correlation is positive, while the 
two-particle one is negative, back-to-back.) With the values of N listed in Table [D] this ranges from 0.2 x 10~ 6 to 
4.8 x 10~ 6 at 158A GeV: for the various centrality bins, this is a factor of 10 smaller than the (vi) 2 V2 values we shall 
find. Therefore, the contamination of correlations due to transverse momentum conservation in our derivation of the 
estimate vi{3} is indeed negligible. 

The flow analysis with either multiparticle method consists of two successive steps. The first step is to estimate 
the average value of V\ and u 2 over phase space (in practice, these are weighted averages, as we shall see shortly), 
which we call "integrated flow". This is done using Eq. JJJ), Eq. © or Eq. 0, which yield v n {2}, w n {4} and «i{3}, 
respectively. The second step is to analyze differential flow, v n (pt , y) , in a narrow (p t , y) window. For this purpose, one 
performs averages as in Eqs. (|7I9|I . where the particle with angle <\>\ belongs to the (pt,y) window under study, while 
the average over cf>2, <j>3, <j>4 is taken over all detected particles. The left-hand sides of Eqs. H7I8JI and the right-hand 
side of Eq. © are then replaced by v n {2}(p t ,y) x v n {2}, v n {4,}{p u y) x w„{4} 3 , vi{S}{p u y) x wi{3}f 2 , respectively. 
This defines the estimates of differential flow from 2, 4 and 3 particle correlations. Note that they can be obtained 
only once the integrated flow is known. In order to reduce the computing time, the analysis was performed over 20 
p t bins of 0.1 GeV/c and 20 y bins of 0.3 rapidity units (instead of 50 y bins in the standard method). 

In the case of NA49, the use of higher order cumulants was limited by statistical errors, in particular for differential 
flow. In practice, up to four estimates of integrated elliptic flow were obtained, namely «2{2}, i>2{4}, t>2{6} and «2{8}, 
but at most two (w 2 {2} and t>2{4}) for differential flow. In the case of directed flow, at most three estimates were 
obtained for integrated flow (wi{2}, i>i{3} and «i{4}) and two for differential flow («i{2} and t>i{3|). 

The practical implementation of these multiparticle methods is described in detail in Refs. |25l l26j. In order to 
illustrate the procedure, we recall here how estimates of integrated (directed or elliptic) flow are obtained from the 
first multiparticle method outlined above. One first defines the generating function 
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(Gn(z)) = (ft ^l + ^M( ze -*^+z*^)^, (10) 

where z is a complex variable, and z* its complex conjugate. The product runs over particles detected in a single 
event, and w n (j) = w n (ptj, yj) is the weight attributed to the j-th particle with azimuthal angle <f>j. Angular brackets 
denote an average over events. A similar generating function for the analysis of directed flow from three-particle 
correlations can be found in Ref. [26|. Weights in Eq. 1|10[> are identical to the weights used in the standard method 
(see Eq. ©), namely y in the center of mass for v\, and pt for v-i- As in the standard method, they are introduced 
in order to reduce the statistical error. 

In Eq. (J2J, we use the same value of M for all events in a given centrality bin. This value has been fixed to 80% of 
the average event multiplicity in the bin. The small fraction of events having multiplicity less than M are rejected. 
For the events having multiplicity greater than M, the M particles required to construct the generating function are 
chosen randomly. Alternatively, one could have chosen for M in Eq. @ the total event multiplicity. We have checked 
on a few examples that results are the same within statistical errors. Note that the value of M is much larger for 
the cumulant method (Table ITTll) than for the standard method (Table ITHl . In the standard method, we have seen 
that cuts were performed in order to minimize the azimuthal asymmetry of the detector, resulting in a lower value 
of M . In the cumulant method, such detector effects are taken into account, as will be explained below, so that cuts 
are not necessary. Furthermore, statistical errors are extremely sensitive to M for higher-order estimates, so that it 
is important to use as many particles as possible. 

The cumulants of 2fc-particle correlations, c„{2fc} are then obtained by expanding in power series the generating 
function of cumulants, C n {z), defined as 



Cn(z) = M^Gniz)) 1 



/M 



1 



+°° , l„|2/c 



H^m^lt^. (id 



fc=0 



The average value of the weight squared, (w 2 ^, has been introduced so that the cumulants c n {2k} are dimensionless. 
In practice, the cumulants are obtained from the generating function using interpolation formulas given in Appendix 
1X1 Finally, each cumulant yields an independent estimate of the integrated flow v n : 

c„{2} = v n {2} 2 , c„{4} = - W „{4} 4 , 
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c n {6} = Av n {6} 6 , c„{8} = -33 V „{8} 8 ... (12) 

A similar procedure holds for differential flow. If the cumulant extracted from the data comes out with the wrong sign 
(for instance, a positive number for c„{4}), one cannot obtain the corresponding flow estimate (w„{4}). As we shall 
see in Sec.E] this does occur, most often for central collisions where the flow is small. There are two reasons for this: 
statistical fluctuations, which may be large for multiparticle cumulants [2^; nonflow correlations, in particular for two- 
particle cumulants, which may be opposite to the correlations due to flow (see the effect of momentum conservation 
on directed flow in Sec IV B|) . 

In the last equation, v n denotes the weighted integrated flow, defined as 

-■ - { n -1. (13) 



It is dimensionless, as it should be. We decided to normalize with \J (w%) rather than (w n ) since the weight w n can be 
negative: for a perfect detector, (u>i) = (y) would vanish! One should note that this integrated flow differs from that 
obtained in the standard method, which integrates the doubly differential flow without weights. Since the average 
values in Eq. (| 1 31) are taken over the whole detector acceptance, the integrated v n is a strongly detector-dependent 
quantity, whose absolute value has little physical significance. It is essentially an intermediate step: as explained 
above, one can analyze differential flow only once integrated flow is known. However, we shall see in Sec. EE] that 
the centrality dependence of v n is meaningful. The magnitude of integrated flow also determines the magnitude of 
statistical errors through the resolution parameter Xn = v n \^M, which is essentially the same quantity as for the 
standard method. Using weights increases \ n roughly by a factor of 1.2. In Table ITTTl are presented the multiplicity 
used in the cumulant method and the corresponding x parameters. 



TABLE III: Listed for the two beam energies and six centralities are: M, multiplicity used in the cumulant method, Eq. 110L 
for the reconstruction of the integrated flow; \ n = v n \/~M, the resolution parameter, where v„ is a (rescaled) weighted flow, see 
Eq. 1131 . The multiplicities are larger for the cumulant method than for the standard method because all particles were used, 
including particles in the backward hemisphere. 



Centrality 


1 


2 


3 


4 


5 


6 


158A GeV 














M 


591 


528 


419 


301 


209 


109 


Xi 




0.27 


0.33 


0.39 


0.41 


0.35 




0.27 


0.42 


0.63 


0.67 


0.62 


0.44 


40A GeV 














M 


318 


257 


185 


120 


80 


42 


Xi 




0.003 


0.19 


0.23 


0.16 


0.26 


X'2 




0.45 


0.43 


0.42 


0.40 


0.23 



Although the formalism may at first sight look complicated, its various features make it the simplest to use in 
practice, for several reasons. First, the several estimates are obtained from a single generating function. Second, 
the generating function automatically involves all possible fc-tuplets of particles in the construction of the /c-particle 
cumulants. Last but not least, the formalism can be used even if the detector does not have perfect azimuthal 
symmetry. In this case, Eqs. ()7I9|I no longer hold. Other terms must be added in order to remove the spurious, 
nonphysical correlations arising from detector inefficiencies, and the number of these terms increases tremendously 
as the order of the cumulant increases. With the generating-function formalism, they are automatically included and 
require no additional work. 

When the azimuthal coverage of the detector is strongly asymmetric, further acceptance corrections must be made, 
which amount to a g lobal multiplicative factor 1 in the relations Eqs. (|12|) between the cumulants c n {k} and the flow 
estimates v n {k} |25l l26l|. In the case of the NA49 acceptance, the corrections are negligible at 158A GeV but can 
become significant at 4QA GeV. We present the range of the corrections factors on the reconstructed flow values in 
Table ITVl 



In general, c n {k} depends not only on v n but also on other harmonics v p with p ^ n. However for most detectors, these interferences 
are negligible and this is indeed the case for the NA49 acceptance. 
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TABLE IV: Acceptance correction factors on the reconstructed values of the flow v n {k}. For integrated flow, the range 
corresponds to the various centrality bins. For differential flow, the correction factors are calculated in centrality bin 4 and the 
range corresponds to the various (pt,y) bins. The larger corrections are for the highest p t or y values. 



158A GeV 


w{2} 


w{3} 


v 2 {2} 


v 2 {4} 


40^1 GeV 


w{2} 


«i{3} 


v 2 {2} 


v 2 {4} 


Integrated 


1.00 


1.01 


1.00 


1.00 


Integrated 


1.03 - 1.04 


1.07 - 1.09 


1.03 - 1.04 


1.01 


Differential 


1.00 - 1.08 


1.01 


1.00 - 1.02 


1.01 - 1.06 


Differential 


1.07 - 1.23 


1.10 - 1.17 


1.07 - 1.29 


1.01 - 1.78 



C. Systematic uncertainties 

Various measures are introduced above in order to quantify azimuthal correlations of particles produced in heavy 
ion collisions. These measures are used in the data analysis with the aim to extract information on directed and 
elliptic flow of primary charged pions and protons emitted in the interaction of two nuclei. The measured correlations 
consist, however, not only of the genuine flow correlations, but also from other physical correlations of primary 
hadrons (nonflow physical correlations) as well as correlations introduced by the imperfectness of the measuring 
methods. The issue of nonflow physical correlations is addressed in the description of the methods in Secs. lIV Al IIV Bl 
and, in particular, in Sec. IV El In the following we discuss various sources of detector induced correlations, as well as 
corrections and cuts used to reduce their influence on the results and the systematic uncertainties. 

The geometrical acceptance of the detector is not uniform in azimuthal angle as seen in Fig. [2] top. This effect is 
corrected for in both methods. As the geometrical acceptance of the detector can be probed to high accuracy by the 
particle yields, the systematic uncertainty caused by non-uniform acceptance is small except as noted below. The 
majority of the events selected by the hardware trigger and off-line event cuts (see Sec. IIII Bfl are Pb+Pb collisions. 
However, there is a small (< 5%) contamination in low multiplicity events from collisions of the Pb beam with the 
material surrounding the target foil. A possible bias caused by this contamination was estimated by varying off-line 
selection cuts on the primary vertex position. No influence on the magnitude of the v\ and v 2 was observed. 

About 90% of tracks selected by the track selection cuts (see Sec. IIII B|l are tracks of primary hadrons coming 
from the main interaction vertex. The remaining fraction of tracks originates predominately from weak decays and 
secondary interactions occurring in the detector material. In order to estimate a possible bias due to this contamination 
the cuts on track distance from the reconstructed event vertex (b Xl b y ) were varied as shown in Fig. 0]with little effect 
on the results. 

The efficiency of track reconstruction and track selection cuts depends on track density in the detector and this 
efficiency is the lowest (~ 80%) for central Pb+Pb collisions at 158^4 GeV at midrapidity. The systematic uncertainty 
due to track losses was estimated by varying track selection cuts (see Sec. IIII B~)l . Additionally, data taken at the two 
magnetic field polarities and during two running periods were analyzed separately and the results compared in Fig. [3] 
The systematic variations are largest in near central collisions for V\ (Avi = 0.005 ) and in peripheral collisions for 
v 2 (Av 2 = 0.01). 

The influence of the particle identification procedures on the flow values was probed by changing the energy loss 
criteria within reasonable limits. The resulting variation of v\ and v 2 are below 0.005. The integrations over p t and 
rapidity involve a weighting procedure on the basis of differential cross sections (see SecEJ. These cross sections are 
available at 40^4 GeV only for pions in central collisions, for other cases they were estimated based on data systematics. 
Variation of the cross section weights within reasonable limits results in variations below 0.001 for v 2 and 0.005 for 
v\. 

The azimuthal coverage of the NA49 TPCs is significantly reduced at AQA GeV as compared to 158^4 GeV (see 
Sec. IIII Cfl . Modifications to the standard analysis method were necessary to reduce the systematic errors. The 
validity of the results from the modified method were scrutinized by applying it to the 158 A GeV data. The results of 
this test are stable and in agreement with those of the standard method for mid-central collisions. In near central and 
peripheral collisions large relative differences are taken as estimate of systematic uncertainties. These are significant, 
if low event multiplicity or low flow values make the event plane determination unreliable. 

The results of our study of systematic uncertainties can be summarized as follows. The systematic error of v 2 (pt, y) 
and vi(pt,y) for pions in mid-central and peripheral collisions is 0.005 and 0.002 respectively. For central collisions 
the numbers increase to 0.01 for both. The error for protons for mid-central and peripheral collisions is 0.005 for v 2 
and 0.01 for v±. For central collisions it increases to about 0.03 for both. At 40^4 GeV these errors could be 50% 
larger than those at 158A GeV. Note that the errors plotted in all figures are statistical ones only. 
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FIG. 4: (Color online) Rapidity dependence of pion elliptic flow for different cuts on the track impact parameter. The sets of 
the cuts were chosen to be 1.0 cm for b x and 0.2 cm for b y for the tight cut. For the open case no cut was used. The standard 
cut is 3.0 cm for b x and 0.5 cm for b y . 
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FIG. 5: (Color online) Centrality dependence of results for different data sets taken with the same beam energy, same centrality 
definition, but two different field polarities and beam intensities. The most recent data set was taken after some changes to the 
detector and analysis procedures. 
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V. RESULTS 

Both methods outlined in Sec. IIVI have been applied at 40A GeV and 158A GeV. The double differential flow 
values v n [p t , y) for each harmonic as obtained from the methods were tabulated as a function of p t , y, and centrality. 
Integration of v n (p t ,y) to obtain v n (p t ) or v n (y) values was done by averaging over the integration variable using 
the cross sections of the particles as weights. The cross section values at 158A GeV had been parametrized and 
were available as a macro. Since no cross sections were available at 40A GeV for non-central collisions, the width 
of the pion Gaussian rapidity distribution and the separation of the two proton rapidity Gaussian distributions were 
scaled down by the ratio of the beam rapidities at 40 and 158A GeV. Since we chose larger y bins in the cumulant 
method than in the standard method, the integration over rapidity was not performed over exactly the same y range. 
More precisely, the upper limit is always smaller for the cumulant method because the results for proton flow did 
not seem to be very stable with respect to integration up to high rapidity values. In the cumulant v(pt) graphs the 
indicated y ranges refer to the cumulant results; the reproduced standard method results in these graphs have the 
y ranges indicated in the preceding standard method plots. The results are presented for three centrality bins (two 
successive bins have been combined, weighted with the known cross sections and the fraction of events in each bin, 
see Table [J and also integrated over the first five centrality bins (weighted with the known cross sections and the 
fraction of the geometric cross section for each bin given in Table Q, whic h we call minimum bias. We present elliptic 
flow (Sec. IV A|l . directed flow fSec. lVB|) . minimum bias results (Sec. IV C^ . centrality dependence (Sec. IV D|l . nonflow 
effects (Sec. IV Ej) . and beam energy dependence fSec. IVF)l . The standard method v\ values have been corrected for 
momentum conservation but the cumulant method «i{2} values have not. 

In the graphs of flow as a function of rapidity the points have been reflected about midrapidity and fitted with 
polynomial curves to guide the eye. Please note that we always use the rapidity in the center of mass, and to calculate 
this the nominal laboratory rapidity of the center of mass was taken to be 2.92 at 158A GeV and 2.24 at 40A GeV. In 
the graphs of flow as a function of p t the smooth curves shown to guide the eye were obtained by fitting to a simple 
hydrodynamic motivated Blast Wave model as described in Ref. [33, 01 but generalized to also describe v\ : 

d(j) b cos(n(j)b)I n (at)K 1 (f3 t )(l + 2s„ cos(n</> 6 )) 

v n {p t ) = — 2^ > (14) 

J* d<t>bI(s{oit)K\{l3t){l + 2s„ cos(n0 b )) 

where the harmonic n can be either one or two, where Iq, I n , and K\ are modified Bessel functions, and where 
oit{4>b) = (Pt/Tf)smh(p((f)b)) and (3t(0b) = (mt/Tf) cosh(jo(<^f,)). The basic assumptions of this model are boost- 
invariant longitudinal expansion and freeze-out at constant temperature Tf on a thin shell, which expands with a 
transverse rapidity exhibiting a first or second harmonic azimuthal modulation given by p(4>b) = Po + Pa cos(n0b). In 
this equation, <f> b — <f> — <&rp is the azimuthal angle (measured with respect to the reaction plane) of the boost of 
the source element on the freeze-out hyper-surface |34| . and po and p a are the mean transverse expansion rapidity 
[vq = tanh(po)] and the amplitude of its azimuthal variation, respectively. The parameters in this model are Tf, the 
temperature, po, the transverse flow rapidity, p ai the azimuthal flow rapidity, and s, the surface emission parameter. 
The Tf parameter was fixed and the s parameter was allowed to be non-zero only for pions at 158^4 GeV. The values 
of the fit parameters themselves are not very meaningful because the flow values derived from two-particle correlations 
contain nonflow effects and the values from many-particle analyses have poor statistics, but the fits do provide the 
curves to guide the eye shown in the pt graphs. The data are clearly not boost invariant, but since we use the blast 
wave model only to fit the pt dependence of the flow, it was felt that a more sophisticated model was not warranted. 
When the fits would not converge the points were just connected, giving rise to the jagged lines in some graphs. 

The error bars shown for the standard method are the standard deviation of the data. For the cumulant method, 
they are calculated analytically following the formulas given in 

Refs^mHi. Tables of the data can be found on the 
web at http://na49info.cern.ch/na49/Archives/Data/FlowWebPage/ 

A. Elliptic flow 

1. 158A GeV 

Results from 158A GeV collisions are displayed in Fig. (from the standard method) and are compared to the 
results from the cumulant method in Figs. and [HI 

We first discuss results for mid-central collisions (Figs. and [SI middle) for the following reason: as usual in flow 
analyses, they have smaller errors than peripheral collisions (Figs. an( l HI top) due to the larger multiplicity, and 
also smaller than central collisions (Figs. and [SJ bottom) due to the larger value of the flow. By errors, we mean 
both the statistical error, shown in the figures, and the uncertainty of the contribution of nonflow correlations, which 
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FIG. 6: (Color online) Elliptic flow obtained from the standard method as a function of rapidity (top) and transverse momentum 
(bottom) for charged pions (left) and protons (right) from 15SA GeV Pb + Pb. Three centrality bins are shown. The open 
points in the top graphs have been reflected about midrapidity. Solid lines are polynomial fits (top) and Blast Wave model fits 
(bottom). 



is not known and not shown in the figures. For mid-central collisions, is positive over all phase space for pions and 
protons. As a function of pt, it rises linearly for pions up to 2 GeV/c (Fig. [7| middle left). For protons (Fig. [7\ middle 
right), the rise is slower at low p t (quadratic rather than linear up to 1 GeV/c), but interestingly, v% reaches the same 
value as for pions at 2 GeV/c. All three methods give compatible results within statistical error bars. As a function 
of rapidity, the pion V2 exhibits the usual bell shape (Fig. |H1 middle left), with a maximum at mid-rapidity. This 
maximum, however, is not very pronounced, and vi remains essentially of the same magnitude, between 2% and 3%, 
over four units of rapidities. For protons (Fig. |H1 middle right), the rapidity dependence is similar. Note that v%{y) is 
slightly larger than for pions, although v-^ijPt) was smaller. This can be explained simply: V2(y) is integrated overp t , 
protons have higher average pt than pions, and i>2 increases with p t . For protons, a small discrepancy appears around 
y = between the estimate from the two-particle cumulants («2{2}) and the standard reaction plane estimate, which 
we do not understand, and consider part of our systematic error. 

For peripheral collisions, i>2 is somewhat larger than for mid-central collisions at low pt (Fig. bottom), but 
comparable at high p t . V2(y) is dominated by the low p t region where the yield is larger, hence it is also larger for 
peripheral than for mid-central collisions (Fig.|Sltop). A small discrepancy can be seen around y — 1 for pions between 
w 2 {2} and the standard u 2 (Fig. |H]top left). 

For central collisions, elliptic flow is much smaller (Fig. EJ). As a consequence, four-particle cumulants could not 
be used, due to large statistical fluctuations. Contrary to more peripheral collisions, «2(y) is larger for pions than for 
protons (Figs. El top, and|SJ bottom), for which it is consistent with zero in the available ra pidi ty range. 

We have also compared these results with those of the earlier analysis published in Ref. [13| and updated on the 
collaboration web page Note that the previous analysis used narrower cuts in p t and y, and had much larger 

statistical errors. Results are compatible for protons, but not for pions. In particular V2(pt) is significantly different at 
low p t , where the increase is now much smoother. The earlier analysis was biased by Bose-Einstein correlations 21], 
whose contribution is reduced in the present analysis thanks to the Wi = p t weights used in Eq. J2Jl , rather than the 
unit weights used in the previous analysis. 
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p, (GeV/c) p, (GeV/c) 

FIG. 7: (Color online) Elliptic flow of charged pions (left) and protons (right) from the cumulant method as a function of 
transverse momentum in 158^4 GeV Pb + Pb collisions. Three centrality bins are shown with results from the standard method, 
from cumulants for two-particle correlations (v2{2}), and from cumulants for four-particle correlations (v2{4}). Solid lines are 
from Blast Wave model fits. 




FIG. 8: (Color online) Same as Fig. [7] as a function of rapidity. The open points have been reflected about midrapidity. Solid 
lines are from polynomial fits. Please note the different vertical scales. 
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FIG. 9: (Color online) Elliptic flow obtained from the standard method as a function of rapidity (top) and transverse momentum 
(bottom) for charged pions (left) and protons (right) from 40 A GeV Pb + Pb. Three centrality bins are shown ("central" 
corresponds to centrality bin 2 only, see Table ID . The open points in the top graphs have been reflected about midrapidity. 
Solid lines are polynomial fits (top) and Blast Wave model fits (bottom). 



2. 40A GeV 

Results from 40 A GeV collisions are displayed in Figs. (from the standard method) . I1UI and llll (from the cumulant 
method). Both the values of V2 and the multiplicities are smaller than at 158A GeV, which results in larger errors. 

In mid-central collisions, v^ipt) (Fig.^|middle) is roughly the same as at 158/1 GeV (Fig. [J| middle), both in shape 
and magnitude, and all three methods give consistent results. The value of v 2 at a fixed p t is in fact very similar up 
to the highest RHIC energy, as pointed out by Snellings 35]. Since the average p t is smaller at AO A GeV than at 
158/1 GeV, however, the corresponding 1*2(2/) is smaller (Fig. ^] middle). As a function of rapidity both the standard 
V2 and «2{2} of pions are remarkably flat. On the other hand, i>2{4} is bell-shaped, as at 158A GeV, but statistical 
error bars prevent any definite conclusion. For protons (Fig. II l| middle right), a discrepancy appears between «2{2}, 
which is bell-shaped, and the standard V2, which has a dip at mid-rapidity. 

For peripheral collisions, V2 is slightly larger, but similar to mid-central collisions (Fig. EJ). The four-particle 
cumulant result i>2{4} could not be obtained due to large statistical errors. 

At this energy, statistics did not allow the determination of the second harmonic event plane, nor the derivation of 
estimates from two- or four-particle cumulants for the most central bin. Thus, the "central" bin corresponds to bin 2 
only in Figs. I§l to llll A striking discrepancy appears between the two-particle estimates (standard V2 and V2{2}) and 
the four-particle cumulant result i>2{4} in the pt dependence (Fig. ^| bottom). It is closely related to that observed 
in the rapidity dependence: For pions (Fig. ^2 bottom left), the standard v 2 and i>2{2} both seem to show a dip at 
midrapidity, where they are compatible with zero; for protons (Fig. ^2 bottom right), the standard V2 and i>2{2} are 
even negative at midrapidity. By contrast, the four-particle cumulant result, t>2{4}, is positive for both pions and 
protons and has the usual bell shape, with a maximum value at midrapidity. 

Since the difference between two-particle estimates and t>2{4} is much beyond statistical error bars, this is a hint 
that the two-particle estimates for pions are affected by nonflow effects (recall that 1*2 {4} is expected to be free from 
two-particle nonflow effects). This could be due to correlations from p deca ys, since p mesons are more concentrated 
near midrapidity, and produce a negative correlation in the second harmonic [22| , which tends to lower the two-particle 
V2 estimates. This effect, however, cannot explain the difference for protons. 

While the four-particle cumulant result «2{4} for central collisions looks very reasonable in shape (bell shape in 
rapidity, regular increase with p t for both pions and protons), its magnitude is unusual: it is as large as for mid-central 
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FIG. 10: (Color online) Elliptic flow of charged pions (left) and protons (right) obtained from the cumulant method as a function 
of transverse momentum from 40^4 GeV Pb + Pb. Three centrality bins are shown. Here, "central" (bottom) corresponds only 
to centrality bin 2 (see Table Pi. Results are shown from the standard method, from cumulants for two-particle correlations 
(«2{2}), and from cumulants for four-particle correlations (i>2{4}). The smooth solid lines are from Blast Wave model fits. 




FIG. 11: (Color online) Same as Fig. 1101 as a function of rapidity. The open points have been reflected about midrapidity. 
Solid lines are from polynomial fits. 



collisions, while we would have expected a significantly smaller value following the observations at 158A GeV. 
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B. Directed flow 

1. 158A GeV 
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FIG. 12: (Color online) Standard directed flow as a function of rapidity (top) and transverse momentum (bottom) for charged 
pions (left) and protons (right) from 158^4 GeV Pb + Pb. Three centrality bins are shown. The open points in the top graphs 
have been reflected about midrapidity. On the top, solid lines are polynomial fits. On the bottom, smooth solid lines are Blast 
Wave model fits. 



Results from collisions at 158 A GeV are shown in Figs. El (from the standard method), El and 1141 (from the 
cumulant method). It can be seen in Fig. El top left that the curves for the different centralities all cross zero 
at midrapidity, indicating that the correction for global momentum conservation in the standard method, shown 
in Fig. for minimum bias data, also works for the individual centralities. One clearly sees the magnitude of the 
correction for momentum conservation in Figs. ll3l and 1 141 by comparing the results from the two-particle cumulants, 
t>i{2} (squares), which are not corrected, with the results of the standard reaction plane analysis (circles), which are 
corrected (see Sec. IIV A"|) . The difference between the two results rises linearly with p t , and is larger for peripheral 
collisions (Fig. El top) than for mid-central collisions (Fig. El middle) , as expected from the discussion in Ref. [23| . 
For central collisions, the negative correlations due to momentum conservation become larger in absolute value than 
the positive correlations due to flow (because flow is much smaller in central collisions), so that the cumulant ci{2} in 
Eq. i|12|) becomes negative and the flow estimate vi{2} could not be obtained. As a function of rapidity, the difference 
between i>i{2} and the standard v\ is approximately constant (Fig. ll4l top and middle). Unlike the standard v\, v\{2} 
does not cross zero at midrapidity, while the true directed flow should. This is a direct indication that «i{2} is biased 
by global momentum conservation (see Fig. [3J top). It is shown merely as an illustration of this effect, and will not 
be discussed further. 

We first discuss the directed flow of pions for mid-central collisions. Its transverse momentum dependence (Fig. 1131 
middle left) is peculiar. The standard v\ (circles) is negative at low pt, but then increases and becomes positive 
above 1.4 GeV/c. The result from the mixed three-particle correlation method [Eq. ©], ui{3}, which is expected to 
be free from all nonflow effects including momentum conservation, is also shown (triangles). Above 0.6 GeV/c, it is 
compatible with the standard reaction plane estimate, but not with the two-particle estimate fi{2}. This suggests 
that nonflow effects are dominated by momentum conservation in this region. Unfortunately, statistical error bars 
on i>i{3} are too large at high p t to confirm the change of sign observed in the standard V\. Below 0.2 GeV/c, 
where momentum conservation is negligible, both the standard analysis value and wi{2} seem to intercept at a finite 
value when p t goes to zero. This is suggestive of nonflow effects arising from quantum Bose-Einstein effects between 
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FIG. 13: (Color online) Directed flow of charged pions (left) and protons (right) from cumulant method, as a function of 
transverse momentum in 158A GeV Pb + Pb. Three centrality bins are shown. Here, "central" (bottom) corresponds to 
centrality bin 2 (see Table U} only. Results are shown from the standard method, from cumulants for two-particle correlations 
(«i{2}), and from cumulants for three-particle correlations (i>i{3}). The ^i{2} estimates are not corrected for momentum 
conservation, while the standard method results are. Smooth solid lines are Blast Wave model fits. Please note the different 
vertical scales. 



identical pions ;21]. The three-particle estimate «i{3}, which is free of nonflow effects, smoothly vanishes at p t = 0. 
This discussion shows that results on directed flow must be interpreted with care, and are more biased by nonflow 
effects than results on elliptic flow. This is mostly due to the smaller value of directed flow (typically 2% in absolute 
value, instead of 3% for elliptic flow). 

The rapidity dependence of the pion v\ for mid-central collisions is displayed in Fig. 1141 middle left. One notes 
that «i{3} vanishes at midrapidity (unlike «i{2}), which confirms that it is automatically corrected for momentum 
conservation effects. Both fi{3} and the standard v\ exhibit a smooth, almost linear rapidity dependence. They 
are in close agreement near mid-rapidity, where we observe clear evidence that the slope is negative. At the more 
forward rapidities, however, the standard v\ becomes larger in absolute value than t>i{3} and seems to saturate, while 
the slope becomes steeper for ui{3}. This small discrepancy can be attributed to the above mentioned Bose-Einstein 
effects which bias the standard v\ at low p t . This bias is only present over the phase space used to determine the 
event plane [21| . i.e. above center-of-mass rapidity y = 1.08. Note that the value of fi{3} is only 4% in absolute 
value at y = 2. This is smaller by at least a factor of 5 than the value obtained by the WA98 collaboration in the 
target fragmentation region (center-of-mass rapidity around -3), where v\ reaches 20% for pions |3f|. There is no 
overlap between their rapidity coverage and ours, so that we cannot check whether the two analyses are consistent. 
Nevertheless, this comparison suggests that the slope becomes much steeper toward beam rapidity, a trend already 
seen on our «i{3} result. 

In peripheral collisions, the pion v\ (Fig. 1131 and Fig. ^] top left) has a p t and y dependence very similar to the 
one in mid-central collisions. Its magnitude, however, is larger, and the increase is more significant than for elliptic 
flow: «i{3} is 50% larger at the most forward rapidities (-6%, instead of -4% for mid-central collisions). In central 
collisions, fi{3} could not be obtained due to large statistical errors. The standard v\ (Fig. 1131 and Fig.^Jtop left) is 
largest at low p t (below 0.3 GeV/c) and forward rapidities (above y = 1) where it is biased by Bose-Einstein effects as 
explained above. This bias is even more important for the centrality bin 1, where a tighter p t cut (below 0.3 GeV/c) 
was chosen for the event plane determination. 

We now discuss the directed flow of protons in mid-central collisions. As a function of pt (Fig.ll3lmiddle right) both 
the standard v\ and t>i{3} are generally positive but almost consistent with zero. The rapidity dependence (Fig. I14l 
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FIG. 14: (Color online) Same as Fig. 1131 as a function of rapidity. The open points have been reflected about midrapidity. 
Solid lines are polynomial fits. 



middle right) is more interesting. v\ is flat near midrapidity (error bars are large, but we can safely state that it is 
flatter than for pions). Only one point, at the most forward rapidity, clearly deviates from zero. This is an essential 
point, since we use it to resolve the overall sign ambiguity of the Vi analysis: at high energies, v\ is assumed to be 
positive for protons at forward rapidities, and this determines the pion flow to be negative in this region (Fig. 1141 
middle left). In peripheral and central collisions where the sign of the proton v\ could not be clearly determined, the 
pion vi was assumed to be negative for sake of consistency. This sign will be established more firmly at 40A GeV (see 
below). As already noted for pions, the slope of v\ above y = 2 must be very steep in order to match WA98 results, 
which give a proton vi of the order of 20% around y = 3 |3|j . 

In peripheral collisions, there is a discrepancy between the standard v\ and «i{3} for protons, which is clearly seen 
on the p t graph (Fig. H3l top right). This difference may be due to nonflow correlations from A decays into protons 
and pions [23, which are automatically corrected for in the three-particle analysis, but not in the reaction plane 
analysis. This contamination from correlations from A decays would then explain why the standard v\ has a negative 
sign at forward rapidities (Fig. ll4l top right). In central collisions, where only the reaction plane estimate is available, 
vi(pt) is compatible with zero (Fig. EHbottom right), but is positive at forward rapidities (Fig. ^] bottom right), as 
expected. 

Let us briefly compare the present results with those of the earlier analysis [l3l Il4| . As in the case of elliptic 
flow, the most significant differences are seen for pions, where v\ was biased in the earlier analysis by Bose-Einstein 
correlations and global momentum conservation |2ll l22j. The earlier vi(p t ) is in fact similar to the present fi{2} 
(squares in Fig.^3 middle left), which suffer from the same biases. 



2. 40A GeV 

Results from 40A GeV collisions are shown in Fig. ^| (from the standard method). 1161 and!17l (from the cumulant 
method). In Fig.^J top, one can see that the standard v\ crosses zero at midrapidity for all centralities, which shows 
that the correlation from global momentum conservation has been properly subtracted. The three-particle results 
t>i{3}, which is automatically free from all nonflow effects including momentum conservation, also crosses zero at 
mid-rapidity for peripheral and mid-central collisions (Fig. 1171 top and middle) . It could not be obtained for central 
collisions due to large statistical fluctuations. The effect of momentum conservation is even larger than at 158^4 GeV, 
as can be seen by comparing v\{2} (not corrected) with the standard v\ (corrected) for peripheral collisions (Figs. El 
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FIG. 15: (Color online) Standard directed flow as a function of rapidity (top) and transverse momentum (bottom) for charged 
pions (left) and protons (right) from iOA GeV Pb + Pb. Three centrality bins are shown. The open points in the top graphs 
have been reflected about midrapidity. On the top, solid lines are polynomial fits. On the bottom, smooth solid lines are Blast 
Wave model fits. 



and 1171 top). For mid-central and central collisions, the momentum conservation effect is so large that wi{2} could 
not even be obtained, as for central collisions at 158^4 GeV. 

The directed flow of pions (Figs. ITBl 1161 and 1171 left) is similar at 40A GeV and at 158 A GeV, both in magnitude 
and shape. The standard v± and «i{3} are compatible. Error bars are larger at the lower energy due to the lower 
multiplicity, especially for the three-particle result «i{3}. The nonzero value of the standard v\ at low p t (Fig. 1X51 
bottom left) probably results from nonflow Bose-Einstein correlations, as at 158^4 GeV. 

The directed flow of protons, on the other hand, (Figs. IT51 1161 and 1171 right), is significantly larger at 40A GeV 
than at 158 A GeV. For mid-central collisions, both the standard V\ and fi{3} clearly differ from zero at forward 
rapidities (Fig. I17| middle right). As already explained, the proton v\ is assumed to be positive at forward rapidities, 
and this fixes the sign of the pion flow to be negative in this region (Fig. El middle left). This is consistent with our 
prescription at 158 A GeV. The standard V\ is larger than «i{3}. The discrepancy is beyond statistical error bars 
(Fig. 1161 middle right) and might be due to some nonflow effect. 

Unlike the pion v±, the proton v± does not seem to be larger for peripheral collisions than for mid-central collisions. 
In Figs. 1121 and 1151 top right panels, the peripheral data seem to exhibit a "wiggle" such that the proton v\ has a 
negative excursion. Due to the large statistical error bars, this can neither be confirmed nor invalidated by the three- 
particle cumulant results in Figs. 1141 and 1171 top right. Nevertheless, such a behavior has been predicted [37| due to 
the variation in stopping in the impact parameter direction in peripheral collisions coupled with the space-momentum 
correlations of flow |38j . This is the first experimental observation of this phenomenon. 



C. Minimum Bias 



The results of the standard method integrated over the first five centrality bins weighted with the fraction of the 
geometric cross section for each bin given in Table|TJare shown in Figs. 1181 andlliJI for the two beam energies. For 40^4 
GeV i>2 bin 1 was not included because we had no results for it. (Also, the cumulant data could not be summed for 
minimum bias graphs because too many centrality bins were missing.) In the lower left graphs, v± for pions shows a 
sharp negative excursion in the first 100 MeV/c of p t . This can also be seen in the graphs for the individual centralities. 
To describe this feature the Blast Wave fits require a very large po parameter. However, the physical explanation is 



23 




FIG. 16: (Color online) Directed flow of charged pions (left) and protons (right) from the cumulant method, as a function of 
transverse momentum in 40A GeV Pb + Pb. Three centrality bins are shown. Results are shown from the standard method, 
from cumulants for two-particle correlations (ni{2}), and from cumulants for three-particle correlations (ui{3}). Please note 
that the «i{2} estimates are not corrected for momentum conservation, while the standard method results are. Smooth solid 
lines are Blast Wave model fits. 
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not clear and the effect may in fact be caused by some very low pt short range nonflow correlation, most probably 
quantum correlations between pions |2lL . Similarly, the positive values of V\ for pions at high transverse momenta 
is due to nonflow correlations (see fi{3} in Fig. Ilfijl. In Fig. 1181 lower right, the proton v\ is consistent with zero at 
all pt values because of the accidental cancellation of its positive and negative values at the different centralitics. 





FIG. 18: (Color online) Standard minimum bias directed and elliptic flow as a function of rapidity (top) and transverse 
momentum (bottom) for charged pions (left) and protons (right) from 158A GeV Pb + Pb. Shown are vi (circles) and vi 
(squares). The open points in the top graphs have been reflected about midrapidity. Solid lines are polynomial fits on the top 
and Blast Wave model fits on the bottom. 



D. Centrality dependence of integrated flow values 

Results from the standard method for the doubly-integrated v n as a function of centrality are shown in Fig. [201 As 
the differential flow values v n (p t ) or v n (y), these results were obtained by averaging the tabulated v n (pt,y) values, 
here over both transverse momentum and rapidity, using the cross sections as weights. For pions on top, the values 
generally increase in absolute magnitude in going from central to peripheral collisions. However, for protons on 
bottom, the values appear to peak at mid-centrality. 

Figure |2~T1 shows the weighted integrated flow values for charged particles, Eq. QliSfl. from the cumulant method 
as a function of centrality, up to eight-particle correlations. Contrary to the standard method results, these values 
were not obtained by averaging the differential values a posteriori. As mentioned in Sec. IIV Bl integrated flow is the 
first outcome of the cumulant method, which is then used to obtain differential values. In particular, this explains 
why for v-i we could derive integrated values from up to eight-particle cumulants, while only four-particle results are 
available for differential flow. Also for V\ we could obtain integrated values from up to four-particle cumulants while 
only three-particle results are available for differential flow with reasonable error bars. 

Please note that the plotted quantity is not v n , but a weig hted flow (tu n e in ^-* JU> )) /y/{w n 2 ). Since it is only 
intended as a reference value in the cumulant method, the weights were chosen so as to maximize it, and indeed it is 
about a factor of 1.2 larger than the standard analysis value because of the weighting. 

A striking feature is the consistency between the flow estimates using more than two particles. Thus, in Fig. 1211 
top left for 158A GeV, t>2{4}, «2{6}, and i>2{8} are in agreement within statistical error bars. This is a clear signal 
that these estimates indeed measure the "true" t>2, and are all arising from a collective motion. In addition, it is 
interesting to note that these three estimates differ significantly from the two-particle estimate i>2{2}, which suggests 
that the latter is contaminated by nonflow effects. This will be further examined in Sec. IV El Apart from that, these 
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FIG. 19: (Color online) Same as FigEHJfor 40 A GeV Pb + Pb. 
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FIG. 20: (Color online) Directed and elliptic flow as a function of centrality for charged pions (top) and protons (bottom) from 
158A GeV Pb + Pb (left) and 40^4 GeV Pb + Pb (right) from the standard analysis. Shown are v± (circles) and «2 (squares). 
The more central collisions are on the left side of each graph. Solid lines are polynomial fits. 



charged particle values follow the usual trend for pions: the absolute value of v\ increases in going to more peripheral 
collisions; i>2 also increases in going from central to semi-central collisions, but then starts to decrease for the most 
peripheral bins, both at 40 and 158A GeV. 
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FIG. 21: (Color online) Weighted flow {w„e 



}/y(wn), Eq. (1131 . with wi = y in the center of mass frame and 



^2 = Pt, from the cumulant method as a function of centrality in 158^1 GeV Pb + Pb (left) and AO A GeV Pb+Pb (right). 
The more central collisions are on the left side of each graph. The lines are polynomial fits. The unplotted points could not be 
obtained or had error bars which were too large. 



E. Nonflow effects 



The purpose of the cumulant method is to remove automatically nonflow effects, so as to isolate how. Nevertheless, 
one may wonder whether this removal is really necessary, and if nonflow effects are significant. 

A first way to estimate the contribution of nonflow effects, which was proposed in Ref. [2(j, consists in plotting 
the quantity N ■ [v n {2] 2 — v n {k} 2 ) as a function of centrality, where N is the mean total number of particles per 
collision in a given centrality bin, and k > 2. The reason is that the two-particle estimate v n {2} is contaminated 
by nonflow effects, while the estimate from more particles v n {k} is not, hence their difference should be due to the 
nonflow correlations. More precisely, the two-particle cumulant reads 

c„{2} = v n {2} 2 

= ( v n) 2 + nonflow 

= K) 2 + |, (15) 

where we have recalled the definition of the two-particle flow estimate v n {2}, see Eq. (1121) . and the nonflow term 
scales as 1/N. Since it is hoped that the multi-particle estimates reflect only flow, v n ~ v n {k}, a straightforward 
rearrangement shows that 

N • K{2} 2 - v n {k} 2 ) = g n (16) 

should be approximately constant. Note that g\ and g 2 were calculated directly from the cumulants, not from the 
values of the flow estimates: 

91 = N ■ ( Cl {2} - Cl {3}/v 2 ) , 

92 = N ■ (c 2 {2} - V^2~{4}) (17) 

(ci{3} denotes the left-hand side of Eq. and the other cumulants have been defined in Sec. IIVB|l . This explains 
why we can display values for almost all centrality bins, while there might be no corresponding flow values because of 
a wrong sign for the cumulant 2 . In Fig. 1221 are shown the coefficients gi and g 2 as a function of centrality, for both 



2 When C2{4} > 0, which occurs for centrality bin 1 at 158^4 GeV, we set C2{4} = in Eq. 1171 . 



27 



158 and 40 A GeV, where the values of N are taken from Table 1TT1 The magnitude of the statistical error follows the 
analytical formulas given in Appendix [H] 
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FIG. 22: (Color online) Nonflow azimuthal correlations from Eq. ll6ll . for the first, g%, (bottom) and second, <?2, (top) Fourier 
harmonics, from 158A GeV (left) and 40/1 GeV (right) Pb+Pb collisions. For gi, the solid points represent all nonflow effects, 
while the open points are corrected for momentum conservation. The horizontal lines are at the mean values. 



As could be expected from the fact that in Fig.|^top left the two-particle points are above the multi-particle points, 
a clear nonflow signal can be observed in Fig.|22ltop left for the second harmonic in full energy collisions. Moreover, 
cj2 is approximately constant, as expected for nonflow correlations: one has approximately i>2{2} 2 — U2{4} 2 + 0.26/iV, 
at least for the four more central bins. 

For the first harmonic, two types of results are displayed in Fig. [23 bottom: the solid points correspond to all 
nonflow effects in Eq. I|ltj|) . while the open points are corrected for momentum conservation. The uncorrected points 
are approximately constant and close to zero at 158,4 GeV for all but the most peripheral bin. This is quite surprising 
since we have seen above that the two-particle fi{2} is strongly contaminated by correlations arising from momentum 
conservation. Actually the contribution of this nonflow effect can be explicitly calculated using Eq. (4) of Ref. [23| : 

mom — cons ill Pt) N 2 q\ 

(y ) (P* )all M 

where averages are calculated over the phase space covered by the detector, except for the subscript "all" which means 
that all particles are taken into account, including those that fall outside the detector acceptance. The corresponding 
quantity (pt 2 ) all is estimated from a model calculation, see Table ITU In the second version of the equation, N is the 
multiplicity over all phase space, M is the multiplicity of detected particles, and / is defined by Eq. (19) of Ref. pjjj . 
Note that for a perfect acceptance, or an acceptance symmetric about midrapidity, ^ mom ^ cons vanishes. As expected, 
^mom-cons depends weakly on centrality and its mean value is about -0.067 at 158A GeV, and -0.10 at 40,4 GeV. 
Since the acceptance is better at 158 A GeV than at 40 A GeV, we expect to get a lower value at that energy, and this 
is the case indeed. Finally, subtracting this contribution to g\ (which amounts to an addition, since ^™ om_cons < o), 
all values are shifted upward and we get a clear positive value at 158 A GeV for all nonflow effects except momentum 
conservation. 

The important point is that the whole gi is only an integrated quantity, averaged over phase space and summed 
over different particle types, and variations of opposite sign may cancel. As a matter of fact, we have already noted 
in Sec. IV Bl that for pions the three-particle estimate fi{3} is lower than i>i{2} at large pt, but larger at low pt 
(Fig. 1131 left): ascribing the large-momentum discrepancy to momentum conservation, and the low-pt one to short- 
range correlations, it turns out that they compensate when one averages over p t . Similarly, we have already mentioned 
that the difference seen in the directed flow of protons (Fig. ED right) between the standard analysis value (corrected 
for p t conservation) and the three-particle estimate at high momentum could be due to correlations from A decays, 
an effect that may explain the open points in Fig l22l bottom. 

On the contrary, the nonflow correlations in the second harmonic, which are seen in Fig. 1221 cannot be localized 
in any definite region of phase space in Figs. 0or|Hl and thus there is no clue regarding the actual effects which 



28 



contribute. One may simply notice that the order of magnitude of 172 is the same as the estimate for the contribution 
from p-meson decays |22| . although resonance decays may only explain part of g 2 . Nevertheless, one can safely 
conclude that nonflow effects are significant at 158A GeV and play a role in the analysis of both V\ and v 2 . 

At 40A GeV, Fig. [23 is less conclusive than at 158^4 GeV, but we have seen in Sec. clear indications of nonflow 
correlations. First, we recall that the importance of nonflow effects, and especially of momentum conservation, in the 
first harmonic is such that it is not possible to obtain directed flow from the two-particle cumulant, «i{2}, while the 
standard method (including correction for p t conservation) or the three-particle cumulant do give results. Apart from 
this large effect, we also mentioned the possible presence of short-range correlations at lowp*, which could explain why 
the standard method value for v\ does not go smoothly to zero at vanishing momentum, while the three-particle result 
Wi{3} does. In the second harmonic, we have seen in Sec. IV Al large discrepancies between the two- and four-particle 
estimates, in central collisions. As suggested above, the observed differences may arise from two-particle correlations 
due to p-meson decays. 



F. Beam energy dependence 
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FIG. 23: (Color online) A comparison of directed and elliptic flow of charged pions at 40 and 158A GeV using the standard 
method is shown as a function of scaled rapidity. The left graph for vi is for peripheral collisions and the right graph for V2 is 
for mid-central collisions. The open points have been reflected about midrapidity. 



A direct comparison between flow at 40 and 158 A GeV is presented in Fig. [53 where v\ and v 2 for pions are plotted 
as a function of scaled rapidity y/ybeam- The use of scaled rapidity is justified by the fact that the width of the 
rapidity distribution of produced hadrons increases approximately as ybeam [ill and that the shapes of the transverse 
momentum distributions are similar at 40 and 158 A GeV Dependence of v\ and v 2 on scaled rapidity is similar 
for both energies, however there is an indication that v 2 is slightly larger at 158A GeV than at AQA GeV, although 
Vi appears to be the same at the two beam energies. For v\ at the other centralities the agreement is not as close, 
but the i>i{3} integrated results, which are free from nonflow effects, agree at the two beam energies. 

Elliptic flow was recently measured at midrapidity in Au+Au at RHIC energies These measurements 

together with the corresponding measurements presented in this paper and other SPS measurements 0> (there 
are newer measurements not taken into account yet, as they were taken at different centralities 0])> as well as the 
measurement at the AGS [48j, allow us to establish the energy dependence of v 2 in a broad energy range. The STAR 
data have been scaled down to take into account the low p t cut-offs of 150 MeV/c at ^snn = 200 GeV and 75 
MeV/c at ^/snn = 130 GeV. The correction factors have been obtained by extrapolating the charged particle yield 
to zero pt and assuming a linear dependence of v 2 (pt) at small transverse momenta. The factors were 1.14 for the 
y/SNN = 200 GeV data and 1.06 for the ^/snn = 130 GeV data. The SPS data do not have a low p t cutoff. In Fig. l2~4l 
v 2 at midrapidity is plotted as a function of center of mass energy per nucleon-nucleon pair for mid-central collisions. 
The rise with beam energy is rather smooth. 

v 2 divided by the initial eccentricity of the overlap region, e, is free from geometric contributions |4?]|. It is useful 
to plot this quantity versus the particle density as estimated by dN/dy of charged particles divided by the area of 
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FIG. 24: (Color online) Energy dependence of V2 near mid-rapidity (0 < y < 0.6 for A0A GeV and < y < 0.8 for 158^4 
GeV) for mid-central collisions (approx. 12 to 34 % of geometrical cross section). The results of NA49 pion vi are compared 
to charged particle v 2 measured by E877, STAR, PHENIX and PHOBOS. 
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FIG. 25: (Color online) v%/e as a function of particle density. The v 2 values are for near mid-rapidity (0 < y < 0.6 for 40A 
GeV and < y < 0.8 for 158^4 GeV). The results of NA49 pion V2 are compared to charged particle vi measured by E877 and 
STAR. The meaning of the horizontal lines (hydro limits) and of the arrow will be discussed in Sec. I VII 



the overlap region, S [38|, |48j . The initial spatial eccentricity is calculated for a Woods-Saxon distribution with a 
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wounded nucleon model from 

_ (y 2 



(19) 



where x and y are coordinates in the plane perpendicular to the beam and x denotes the in-plane direction. Fig. 
shows our results for t>2{2} together with the recent STAR results on elliptic flow from 4-particle cumulants |39ll4£ 
and E877 results from the standard method. The RHIC data have been corrected for their low p t cut-off as 
described for Fig. [21 Only statistical errors are shown. 

The NA49 results shown in the figure were obtained with the 2-particle cumulant method. The total systematic 
uncertainties for the points are unfortunately rather large and amount to about 30% of the presented values. The 
systematic errors for the most central collisions are even larger and taking into account also the larger statistical 
errors we do not plot them in the figure. The results obtained for NA49 with the standard method as well as from 
4-particle cumulants have larger statistical errors but are within the systematic errors quoted. There are also possible 
systematic errors of the order of 10-20% due to uncertainties in the centrality measurements from one experiment to 
the other. Note that the most peripheral RHIC data are slightly above the slightly more central data from the SPS, 
but would agree within these systematic errors. This slight non-scaling could be related to the use of e calculated 
with a weight proportional to the number of wounded nucleons. A different weight could give somewhat different 
relative values for spatial eccentricities of central and peripheral collisions. In general the figure shows a nice scaling 
of the elliptic flow divided by the initial spatial eccentricity when plotted against the produced particle density in the 
transverse plane. The physical interpretation will be discussed in Sec. IVII 

The energy dependence of directed flow is also instructive. Directed flow has not yet been seen at RHIC, but it has 
been extensively studied at lower energies at the AGS. The most striking difference between the present results and 
AGS results is in the rapidity dependence of the proton directed flow. Up to the top AGS energy, the proton v\{y) 
follows the well-known S-shape, with a maximum slope at midrapidity |4^. By contrast, at SPS, already at AO A GeV, 
the slope at midrapidity is consistent with zero within our errors. This is compatible with a smooth extrapolation of 
AGS results, which already show a significant decrease of the slope at midrapidity between 2 and 11 A GeV. 

While the proton v\ near midrapidity becomes much smaller as the energy increases, the pion v\ near midrapidity 
remains of comparable magnitude. The main difference, compared to AGS energies, is that v\ remains negative until 
very high values of pt, while at AGS it becomes positive typically above 500 MeV/c ^fji which was interpreted as an 
effect of the sidewards motion of the source. 



VI. MODEL COMPARISONS 



Here we review theoretical predictions for elliptic and directed flow at SPS energies and compare them with our 
results. Note however that detailed predictions have only been made for 158A GeV collisions. 



A. Elliptic flow 



Elliptic flow at ultrarelativistic energies is interpreted as an effect of pressure in the interaction region. In the 
transverse plane, particles are created where the two incoming nuclei overlap. This defines a lens-shaped region 
for non-central collisions. Subsequent interactions between the particles drive collective motion along the pressure 
gradient, which is larger parallel to the smallest dimension of the lens. This creates in-plane, positive elliptic flow Q. 
At early times, however, spectator nucleons tend to produce negative elliptic flow |52j : while this effect dominates at 
energies below 5 GeV per nucleon |5^|, one expects it to be negligible at SPS energies, except close to the projectile 
rapidity, which is not covered by the present experiment. 

As it is essentially a pressure effect, elliptic flow is sensitive both to the equation of state of nuclear matter (i.e., the 
relation between pressure and density), and the degree of thcrmalization reached in the system. The equation of state 
enters mainly through the velocity of sound, c s = (dP/de) 1 / 2 , which controls the magnitude of pressure gradients: 
naturally, a "softer" equation of state, with smaller c s , produces smaller elliptic flow, and this is true in particular 
in the presence of a first-order phase transition where c s vanishes From a microscopic point of view, elliptic 
flow increases with the number of collisions per particle: very generally, one expects a linear increase, followed by a 
saturation when the system reaches thermal equilibrium [54j. Disentangling both effects (equation of state vs. degree 
of thcrmalization), however, is far from obvious. 

Various predictions have been made for elliptic flow at SPS energies. They are based either on hydrodynamic 
models, where full thermal equilibrium is assumed, or on transport models, where collisions between particles are 
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treated at a microscopic level. We briefly review these predictions and discuss how they compare with the present 
data. 

Let us start with the transverse momentum dependence of V2- As can be seen in Figs. 171 and 1101 V2 vanishes at 
low p t and then increases regularly. The dependence is almost linear for pions, but rather quadratic below 1 GeV/c 
for protons: therefore i>2 is smaller for protons than for pions at low p t , while they become of comparable magnitude 
at 2 GeV/c. These nontrivial features are qualitatively reproduced by hydrodynamical models 0, where the 
dependence on the hadron mass is well understood |34j | . Unfortunately, there do not seem to be any predictions from 
transport models for v^ipt) at SPS ener gies . Several such calculations are available at the higher RHIC energies, 
where the same features are observed [5(| |52I Ha, |53 • Coming to quantitative comparisons, elliptic flow at 158^4 GeV 
seems to be overestimated by hydrodynamical calculations, although not by a large factor. Kolb et al. Q predict 
V2(pt — 1 GeV/c) ~ 12% for pions in mid-central collisions, while our various estimates (see Fig. middle left) are 
closer to 8%. The model of Teaney et al. [6(J, which couples hydrodynamics to a transport model, naturally yields 
smaller values of v-i (about 10% for pions at pt = 1 GeV/c), but still higher than our data. These hydrodynamical 
calculations use "soft" equations of state, which are needed in order to reproduce the measured pt spectra. The fact 
that they tend to overestimate V2 suggests that thermalization is only partially achieved even at the top SPS energy. 

The rapidity dependence of vi is also instructive. Due to the strong Lorentz contraction at ultrarelativistic energies, 
all the produced particles essentially originate from the same point z = t = Q. Therefore, the longitudinal velocity 
v z of a particle (or its rapidity) is expected to be strongly correlated to its position z throughout the evolution of 
the system. This means that particles with the same rapidity are those which interacted throughout the evolution. 
Following the common interpretation of elliptic flow as due to secondary collisions, one expects an increase of i>2 with 
the multiplicity density dN/dy, up to the limiting value corresponding to thermal equilibrium. It is interesting to note 
that, at 158^4 GeV, the rapidity dependence of V2 (Fig. [SJ follows that of the multiplicity dN/dy j32, although it is 
less pronounced. This again suggests only partial thermalization. Let us now come to quantitative comparisons with 
the various models. Most of the calculations were done with the same pt cuts as in the earlier NA49 analysis i.e. 
50 < pt < 350 MeV/c for pions, which would yield values of 1*2(2/) smaller than in the present analysis by a factor 
of approximately 2 according to our measurements. First, all above mentioned hydrodynamical calculations assume 
rapidity independence [6l|, so that they are unable to predict the rapidity dependence of 1)2- Hirano [62( performs a 
full three-dimensional hydro calculation, tuned to reproduce both rapidity and transverse momentum spectra. The 
resulting V2 (y) is about 3% at midrapidity, significantly larger than our result in Fig. [H] middle left (with due attention 
to the different cuts in p t , as explained above), and the decrease at higher rapidities seems to be steeper than in our 
data. Various predictions are available from microscopic models for V2(y) of pions and protons at the top SPS energy. 
Predictions of RQMD |6^| and UrQMD [64| seem to be in agreement with our data, at least qualitatively. Calculations 
based on the QGSM model |6{| show a dip near midrapidity for peripheral collisions, which we do not see in the data. 

The centrality dependence of 1*2 is particularly interesting. In hydrodynamics, it closely follows the centrality 
dependence of the initial asymmetry, i.e., V2 is proportional to the eccentricity e of the lens-shaped object It thus 
increases with impact parameter up to very peripheral collisions if the system is thermalized. On the other hand, 
if only partial thermalization is achieved, i>2 is smaller for peripheral collisions, and the maximum of «2 is shifted 
to less peripheral collisions |iM l55j. Since the eccentricity of the lens-shaped object, e, increases linearly with the 
number of spectator nucleons one would expect a linear increase of t>2 with E° in Figs.l2(Jlandl21l At 158A GeV, 
the experimental curve clearly bends over in the two most peripheral bins, and this tendency is more pronounced 
with estimates from higher order cumulants which are less sensitive to nonflow correlations. This again suggests that 
thermalization is not fully achieved. At 40^4 GeV, error bars are larger, but the decrease of V2 is even more pronounced 
in the most peripheral bin: departure from thermalization would then be more important than at the higher energy, 
which is to be expected since the density is lower. Finally, a peculiar centrality dependence of V2 was proposed as a 
signature of the phase transition to the quark-gluon plasma |47| : with increasing centrality (i.e. decreasing impact 
parameter), the ratio V2/e was predicted to first increase, due to better thermalization, then saturate, and eventually 
increase again for the most central collisions as a result of the re-hardening of the equation of state when a quark-gluon 
plasma is produced. At 158^4 GeV, however, Fig. 1211 shows that the increase of V2 with the number of spectators 
is linear for the most central collisions, which suggests a constant t^/e- At 40 A GeV, on the other hand, the flow 
obtained from 4-particle cumulants remains unusually high for central collisions, which may be a hint of interesting 
physics. 

We finally come to the dependence of V2 on beam energy. There are very few model calculations of V2 at 40A GeV 
in the literature. The most natural expectation is that V2 should increase regularly as the beam energy increases, and 
eventually saturate: first, because the effect of spectators (which tends to lower v 2 ) decreases, second, because the 
degree of thermalization increases. Such a behavior is indeed observed with a smoothly varying equation of state • 
If a strong first order phase transition occurs in the system, however, an interesting structure is predicted 0, where 
a small maximum in V2 occurs somewhere between the top AGS and the top SPS energies. Such a structure does not 
appear in Fig. [2U which shows that the increase of V2 from AGS to RHIC is smooth. 
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There is a nice way to combine the centrality dependence and the beam energy dependence of v 2 , by plotting the 
scaled anisotropy v 2 /e versus the particle density as estimated by dN/dy of charged particles divided by the area 
of the overlap region, S j^S H^- If the system is locally thermalized, and if there is no marked structure in the 
equation of state (such as a strong first order phase transition), v%/e is expected to be centrality independent as 
explained above, and to depend only weakly on the beam energy through the equation of state: this is indicated 
by the horizontal lines ( "hydro limits" ) in Fig. EH If the mean free path of the particles becomes as large as the 
transverse size, thermalization is not reached, and v%je is smaller than the hydro limit: it is expected to increase 
smoothly with the ratio between the transverse size and the mean free path, and eventually saturate at the hydro 
limit. The latter ratio scales itself like the charged particle density per unit transverse area (l/S)dN/dy, which is the 
abscissa in Fig. [23] The data shown in this figure are compatible with the above picture, since data taken at different 
beam energies fall on the same universal curve. It is quite remarkable, and far from obvious, that this simple physical 
picture allows us to compare quantitatively a peripheral collision at the full RHIC energy with a central collision at 
the lower SPS energy. 

A further interest in this plot is the speculation 0] whether there is flattening at an abscissa value of about 10 
which could indicate a change in the physics of rescattering and suggest a deconfinement phase transition in this 
region. It is noteworthy that the color percolation point |50j . which is shown by an arrow in the figure, is just in the 
same region. Unfortunately, taking into account the large systematic errors of the results, we are not able to resolve 
this question. Further data analysis at 80A GeV SPS and v / s]v7V = 20GeV RHIC beam energies would be useful. 

B. Directed flow 

Elliptic flow at SPS is created by the interactions between the produced particles, but is to a large extent independent 
of the mechanism of particle production. By contrast, the physics of directed flow probes a time scale set by the 
crossing time of the Lorentz-contracted nuclei |66| . It is therefore more subtle, and more strongly model dependent 
than the physics of elliptic flow. 

Directed flow defines the impact parameter as an oriented arrow in the transverse plane. As such, it reflects the 
asymmetry between target and projectile (to which elliptic flow is insensitive). For this reason, most model calculations 
only address the directed flow of nucleons. In addition, most papers do not deal with v\, but with the mean transverse 
momentum projected on the reaction plane, (p x ) = (ptvi(pt,y)}, as a function of rapidity, for historical reasons llSl. 

It was first expected that directed flow would be negligible at SPS energies. In 1991, it was predicted et al. 
that it might however be large enough to be measurable in Pb-Pb collisions. Non-zero v\ for nucleons was predicted 
both by a transport model (QGSM) including re-scatterings, and by a hydrodynamical model. Furthermore, in the 
hydro model, (p x ) depended strongly on the equation of state: as elliptic flow, directed flow is smaller with a softer 
equation of state (involving for instance a phase transition to a quark-gluon plasma). 

The interest in directed flow was revived following the prediction that the "softest point" of the equation of state 
could be directly observed at the AGS |6Sl |. A deep minimum of (p x ) may appear, at an energy of about 6 A 
GeV [63 • These predictions, however, crucially rely on the assumption that the early stag es of the collision, when 
particles are produced, can be described by one-fluid hydrodynamics. A two-fluid model [7(j predicts no minimum as 
energy increases. In a three- fluid model 71], on the other hand, a minimum occurs, but at a higher energy, around 
10-20 GeV per nucleon. It is followed by an increase up to a maximum at 40^4 GeV. Unfortunately, no quantitative 
estimate is provided. A similar structure is predicted in the transport model UrQMD |66|. where the minimum of (p x ) 
also appears as a consequence of the softening of the equation of state (although no quark-gluon plasma is explicitly 
incorporated in the model): (p x ) increases from the top AGS energy and then saturates above 40 GeV per nucleon at 
a value of 30 MeV/c. Unfortunately, an observable such as (p x ) is largest in the target and projectile rapidity regions, 
which are not covered by the present experiment, so that quantitative comparisons are not possible. 

We therefore concentrate on predictions for differential flow as a function of rapidity, Vi (y) . The remarkable feature 
of our data in Figs. H4l and ll7l is the flatness of v\ for protons near mid-rapidity, where the slope at midrapidity may 
even be negative for peripheral collisions. This "antiflow," or "wi ggle " was predicted at RHIC energies I37J. It is also 
present at SPS energies in most fluid-dynamical calculations [7(J, UM> UM ■ Transport models like UrQMD |64j predict 
a too large value of v±(y) for protons, although they do see some flattening near midrapidity. In the QGSM model 
|65|. v\(y) for protons shows no flattening, except for peripheral collisions. The value of v\ of protons is much larger 
than our data: about 7% for central collisions and 10% for mid-central collisions at y = 1.5, while we see at most 2%. 
The value of Ui(pt) in QGSM is also larger than our data for protons. 

There are comparatively very few studies of pion directed flow at SPS energies. Pion directed flow was first seen 
in asymmetric collisions |74| at lower energies. The v\ of pions usuall y h as a sign opposite to the proton vi, which is 
understood as an effect of shadowing and absorption on nucleons |75ll76|. The same effect was observed at AGS [?3 
and is clearly seen in all our data. UrQMD [64( correctly predicts vi{y) for pions at the top SPS energy, although 
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their prediction for protons is too high. In the QGSM model [65|, the pion v\{y) has a minimum at about half the 
projectile rapidity. In Fig. ^] the same tendency is seen in the standard analysis (also quantitative agreement) but 
not in the three-particle estimate which shows further decrease toward projectile rapidities. The peculiar transverse 
momentum dependence of v± for pions, which is negative at lowp t and then positive above 1 GeV/c (see Figs. Upl and 
HB)l . does not seem to have been predicted. 

VII. SUMMARY 

We have presented a true multiparticle analysis of directed and elliptic flow which provides results on v\ and V2 
values of pions and protons as function of rapidity, transverse momentum, collision centrality and beam energy (40 
and 158 A GeV). 

Two independent analyses were carried out using two different methods: the standard method of correlating particles 
with an event plane, and the cumulant method of studying genuine multiparticle correlations. The cumulant method 
yields several independent estimates of the flow. From the point of view of physics, the lowest-order estimates 
«i{2}, i>2{2} are essentially equivalent to the estimates from the event-plane method, while higher-order estimates 
i>i{3}, i>2{4} are hopefully free from nonflow effects. This is the main motivation of the cumulant method. The 
two methods are very different in their practical implementation. The cumulant method no longer requires one to 
construct subevents, or to correct for the event-plane resolution. All flow estimates are derived from a single generating 
function of azimuthal correlations. Constructing this generating function, however, requires more computer time than 
the standard flow analysis. Another significant difference between the two methods is that the cumulant method 
takes naturally into account azimuthal asymmetries in the detector acceptance. Hence the flattening procedures, 
and the cuts in phase space which are required in the event-plane method in order to minimize the effects of these 
asymmetries, are no longer required. The price to pay for all these enhancements is increased statistical errors. 

We have obtained the first direct, quantitative evidence for collective motion at these energies: elliptic flow at 158A 
GeV has been reconstructed independently from genuine 4, 6 and 8-particle correlations, and all three results agree 
within statistical errors (Fig. top left). This is confirmed at both energies by differential analyses of elliptic flow 
(as a function of rapidity or transverse momentum) from genuine 4-particle correlations. In the case of directed flow, 
nonflow correlations due to momentum conservation, which are large, have been subtracted. Furthermore, a new 
method of analysis from 3-particle correlations, which is unbiased by nonflow correlations, has been implemented for 
the first time at both energies. 

The directed flow of protons reveals a structure which is characteristic of ultrarelativistic energies, and is not present 
at AGS energies. A clear separation appears for the first time between the central rapidity region, where the proton 
v\ is essentially zero, and the target-projectile fragmentation region, where it is large. Indeed, at AOA GeV, significant 
directed flow is observed only at the most forward rapidities covered by the detector acceptance (Fig. El right). 
At 158 A GeV, where the acceptance covers smaller values of the scaled rapidity, v\ values are consistent with zero 
fFig. 1141 right), within statistical errors and possible contributions by nonflow effects. In the fragmentation region, 
on the other hand, large v\ values have been observed by WA98 |5l|. At both energies, the first observation of the 
"wiggle" (i.e., a negative slope of the proton v\ near mid-rapidity) is reported, but there are indications that it may 
be due to nonflow effects. 

Surprisingly, the directed flow of pions does not follow the same behavior as that of protons. While the proton v\ 
at central rapidity is much smaller than at AGS energies, the pion v\ remains essentially of the same magnitude. It 
becomes even larger, in absolute value, than the proton v\. This amazing phenomenon, which has never been observed 
at lower energies, clearly indicates that the proton v\ and the pion v\ have different physical origins. The directed flow 
of pions behaves similarly at the two beam energies, both in magnitude and in shape. It has a peculiar, essentially 
flat, transverse momentum dependence fFigs. IT31 and 1161 left). Its centrality dependence is also quite remarkable: it 
increases in magnitude steadily without saturating up to the most peripheral collisions. (Fig. 1201 top, and Fig. 1211 
bottom). 

Elliptic flow becomes the dominant azimuthal anisotropy at ultrarelativistic energies. While it is smaller than 
directed flow up to the top AGS energy, here it becomes larger already at 4QA GeV. This is again an indication that 
SPS is probing the truly ultrarelativistic regime. As a consequence of the larger value, our estimates of i>2 are more 
accurate than our estimates of v\. As a function of transverse momentum, V2 increases almost linearly for pions, and 
more quadratically for protons, as already seen at RHIC. On the other hand, the rapidity dependences are the same 
for pions and protons: t>2 is approximately constant in the central rapidity region (see e.g., Fig. El middle) and drops 
in the target-projectile fragmentation regions, roughly where the proton v\ starts increasing. At 158A GeV, v<i has 
a pronounced centrality dependence (Figs. l2"Ulandl?Tl left) but, unlike directed flow, it saturates for very peripheral 
collisions. The centrality dependence at 4(M GeV (Figs. 1201 and [2"T1 right) is less significant. For protons (shown in 
the same figures) the tendencies are the same, but the significance is reduced due to larger statistical errors. 



34 



The energy dependence of v 2 looks roughly linear from AGS up to RHIC energies (Fig. I24fl . No indication of 
non-monotonic behavior is visible, as would be expected from the softening of the equation of state for a system close 
to the critical temperature. The dependence of elliptic flow (divided by the eccentricity of the nuclear overlap region) 
on particle (rapidity) density also exhibits a smooth increase without significant structure (Fig. I25fl which would 
indicate a change in the physics of rescattering. 
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APPENDIX A: INTERPOLATION FORMULAS FOR THE CUMULANTS 



In this appendix, we give interpolation formulas which were used to extract the cumulants of multiparticle corre- 
lations up to eight-particle correlations. Following the procedure proposed in Ref. |25| . we introduce interpolation 
points Zp. q = x p , q + iy p , q , 



cos 



2qir 

~9~ 



Vp,q = r oVP sin ( ~~7T~ ) . ( A1 ) 



for p = 1, • • • , 4, q = 0, • • • , 8, and ro of order unity. The generating function of cumulants Eq. IjllH is then computed 
at these various points: C p . q = C n (z p , q ), and averaged 

1 8 

C p = - ^2 C P ,q- (A2) 
y q =o 

Finally, the cumulants for correlations of 2, 4, 6, and 8 particles are given by 
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c «{ 6 } = R ,\ ,3 t9Ci-12C2 + 7C 3 -|c4), 
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"' !' s ! = TT^ (-Ci + \C 2 - C a + i(7 4 ) . ( A3 1 
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APPENDIX B: STATISTICAL ERRORS OF THE CONTRIBUTION OF NONFLOW EFFECTS 

In this appendix, we derive the statistical uncertainties on the coefficients a, in Eqs. (|1 T|) . using the analytical 
formulas for the statistical errors of the various estimates v n {k} given in Refs. |25ll26l|. 

Let us start with that measures nonflow correlations in the second harmonic. We recall its definition: 

g 2 = N ■ (c 2 {2} - v 7 -^}) , 

= N ■ (v 2 {2} 2 - v 2 {4} 2 ) , (Bl) 
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where the second line holds provided C2{4} has the right, negative sign. The statistical error of g 2 can easily be 
calculated and reads 

(Sg 2 ) 2 = AN 2 v 2 [{8v 2 {2}) 2 + {5v 2 {A}) 2 - 2 ({v 2 {2}v 2 {4}) - (v 2 {2}) («a{4}»] . (B2) 

where we used S(v 2 ) = 2v n 5v n . Note that since v 2 {2} and v 2 {4} are slightly correlated if they measure an existing 
v 2 , the error is smaller than the mere geometrical mean of the uncertainties on each flow estimate. Using formulas 
given in Ref. (2^, we obtain 



6 » = m7n^ ' (B3) 



where N is the total number of emitted particles, M the multiplicity used in the analysis, -/V ov ts the total number of 
events and \2 — v 2 \f~M the resolution parameter. The estimate of v 2 from 4-particle cumulants, i>2{4}, was used to 
calculate \2- 

Eq. 1|B3|1 is valid as long as relative errors on both 2-particle and 4-particle cumulants are small. When relative 
errors become large (usually for central collisions where v 2 is small, or for the most peripheral collisions where M is 
small), it may happen that the cumulant C2{4} is positive, so that g 2 in Eq. ijBlfl is undefined. In this case, we set 
our 4-particle estimate «2{4} to zero in the second line of Eq. (|B1(I . We then estimate the statistical error noting that 
it is essentially dominated by that on the 4-particle cumulant C2{4}, which is <5c2{4} = 2/ (M 2 ^/iV cvts ) 25] when the 
resolution parameter \ 2 is small. The statistical error of g 2 then reads 

NV2 



Sg 2 - N^iA} = v — . (B4) 



Consider now gi, that is, nonfiow correlations in the first harmonic. We may assume that ci{2} and ci{3} are 
uncorrelated since the latter comes from a mixed correlation. Therefore, the statistical error of gi can be written as 



(Sgif = A' 2 



t>2 V V 2 



(B5) 



If the estimated «i{3} is reconstructed from «2{4}, as is the case in the 158A GeV analysis, one should use the 
statistical error (5^2{4} in Eq. (JBBJ. However, at 40A GeV, since we used ^{2} (see Sec. HVB() we have to use 5v 2 \2\ 
in Eq. (|B5|I . The corresponding formulas are given in [23, while the formula for <5ci{3} can be found in Ref. |2q|. 



S. Voloshin and Y. Zhang, Z. Phys. C 70, 665 (1996) arXiv:hep-ph/9407282 . 
W. Reisdorf and H. G. Ritter, Ann. Rev. Nucl. Part. Sci. 47, 663 (1997). 

J.-Y. Ollitrault, Phys. Rev. D 46, 229 (1992). 

D. H. Rischke, Nucl. Phys. A 610, 88c (1996) arXiv:nucl-th/9608024 . 
J.-Y. Ollitrault, Nucl. Phys. A 638, 195c (1998) [arXivinucl-ex /9802005 . 

P. Danielewicz, Nucl. Phys. A 661, 82c (1999) |arXiv:nucl-th/9907098r^ 

P. F. Kolb, J. Sollfrank, and U. Heinz, Phys. Lett. B 459, 667 (1999) arXiv:nucl-th/9906003 . 

~ i!29l- 



P. F. Kolb, J. Sollfrank, and U. Heinz, Phys. Rev. C 62, 054909 (2000) arXiv:hep-ph/0006129 
D. Teaney, J. Lauret, and E. V. Shuryak, Phys. Rev. Lett. 86, 4783 (2001) |arXiv:nucl-th/OOljo58j . 
J. Bachler et al. [NA49 Collaboration], Searching for QCD Phase Transition, Addendum-1 to Proposal SPSLC/P264, 
CERN/SPSC 97 (1997). 

S. V. Afanasiev et al. [NA49 Collaboration], Phys. Rev. C 66, 054902 (2002) arXiv:nucl-ex/0205002 . 
M. Gazdzicki and M. I. Gorenstein, Acta Phys. Polon. B 30, 2705 (1999) |arXiv:hep-ph/9803462f~ 

H. Appelshauser et al. [NA49 Collaboration], Phys. Rev. Lett. 80, 4136 (1998) arXiv:nucl- ex/9711001| ; A. M. Poskanzer, 
Nucl. Phys. A 638, 463c (1998). 

http://na49info.cern.ch/na49/Archives/Images/Publications/Phys.Rev.Lett.80:4136-4140 1998/ . 

A. M. Poskanzer et al. [NA49 Collaboration], Nucl. Phys. A 661, 341c (1999) |arXiv:nucl'ex/9906013|. 

A. Wetzler et al. [NA49 Collaboration], Quark Matter '02, Nucl. Phys. A 715, 583c (2003) |arXiv:nucl-~ /0212023|. 

N. Borghini, P.M. Dinh, and J.-Y. Ollitrault, Quark Matter '02, Nucl. Phys. A 715, 629c (2003) |arXiv:nucl-'tb' /0208014 . 
P. Danielewicz a nd G. Odyniec, Phys. L ett. B 157, 146 (1985). 
J.-Y. Ollitrault, [arXiv:nucl-ex/971 1003] . 

A. M. Poskanzer and S. A. Voloshin, Phys. Rev. C 58, 1671 (1998) arXiv:nucl-ex/9805001 . 



3G 



P. M. Dinh, N. Borghini, and J.-Y. Ollitrault, Phys. Lett. B 477, 51 (2000); 
N. Borghini, P. M. Dinh, and J.-Y. Ollitrault, Phys. Rev. C 62, 034902 (2000). 

N. Borghini, P. M. Dinh, J.-Y. Ollitrault, A.M. Poskanzer, and S.A. Voloshin, Phys. Rev. C 66, 014901 (2002) 

arXiv:nucl-th/0202013 . 

N. Borghini, P. M. Dinh and J. Y. Ollitrault, Phys. Rev. C 63, 054906 (2001 ) |arXiv:nucl-th/0007063| . 

N. Borghini, P. M. Dinh, and J.-Y. Ollitrault, Phys. Rev. C 64, 054901 (2001) arXiv:nucl-th/0105040 ; see also Proceedings 

of the International Workshop on the Physics of the Quark-Gluon Plasma, Ecole Polytechnique, Palaiseau, France, p. 237 

(2001) arXiv:nucl-ex/0110016 . 

N. Borghini, P. M. Dinh, and J.-Y. Ollitrault, Phys. Rev. C 66, 014905 (2002), |arXiv:nucl-t*h /0204017 . 
S. Afanasev et al. [NA49 Collaboration], Nucl. Instrum. Meth. A 430, 210 (1999). 

H. Appelshauser et al. [NA49 Collaboration], Eur. Phys. J. A2, 383 (1998). 

K. H. Ackermann et al. [STAR Collaboration], Phys. Rev. Lett. 86, 40 2 (2001) |arXiv:nucl-ex/00 0901 1 . 
J.-Y. Ollitra ult, Phys. Rev D 48, 1 132 (1993) [arXiv:hep-ph/9303247| . 
N. Borghini, arXiv:hep-ph/0302139 

G. E. Cooper, PhD thesis, University of California at Berkeley, 2000 [LBNL-45467, 
http://na49info.cern.ch/cgi-bin/wwwd-util/NA49/NOTE7222 : G. E. Cooper [NA49 Collaboration]. Nucl. Phys. 
A661, 362c (1999). 

C. Adler et al. [STAR Collaboration], Phys. Rev. Lett. 87, 182301 (2001) [arXiv:nucl-e x /0107003] . 

P. Huovinen, P. F. K olb, U. W. Heinz, P. V. Ruuskanen, and S. A. Voloshin, Phys. Lett. B 503, 58 (2001) 

|arXiv:hep-ph/0101136| . 

R. Snellings [the STAR Collaboration] , |arXiv:nucl-ex /030500T1 

M. M. Aggarwal et al. [WA98 Collaboration], arXiv:nucl-ex/9807004 . 

R. J. Snellings, H. Sorge, S. A. Voloshin, F. Q. Wang, and N. Xu, Ph ys. Rev. Lett. 84, 2803 ( 2000) [arXiv:nud-ex/ 9908001 . 
S. A. Voloshin, Quark Matter '02, Nucl. Phys. A 715, 379 c (2003) |arXiv:nucl -ex/0210014 . 

L. Ray, Quark Matter '02, Nucl. Phys. A 715, 45c (2003) arXiv:nucl-ex/0211030 . 

S. Esumi et al. [PHENIX Collaboration], Quark Matter '02, Nucl. Phys. A 715, 599c (2003) |arXiv:nucl-ex/0210012l. 
S. Manly et al. [PHOBOS Collaboration], Quark Matter '02, Nucl. Phys. A 715, 611c (2003) arXiv:nucl-ex/0210036 . 
K. Filimonov et al. [CERES Collaboration], International Nuclear Physics Conference, Berkeley, (2001) [arXiv: 

nucl-ex/0109017 and private communication. 

J. Slivova [CERES Collaboration], Quark Matter '02, Nucl. Phys. A 715, 615c (2003) arXiv:nucl-ex/0212013 . 
S.I.Esumi, J. Slivova, J. Milosevic for CERES Collaboration SFIN, year XV, Series A: Confe rences, No. A2(2002) 

H. Liu et al. [E895 Collaboration], Phys. Rev. Lett. 84, 5488 (2000) [arXiv:nucl-ex/0005005l . 

J. Barrette et al. [E877 Collaboration], collisions at 11-A-GeV/c," Phys. Rev. C 56, 3254 (1997) arXiv:nucl-ex/9707002 . 

H. Sorge, Phys. Rev. Lett. 82, 2048 (1999) |arXiv:nucl-th/9812057| 

S. A. Voloshin and A. M. Poskanzer, Phys. Lett. B 474, 27 (2000) arXiv:nucl-th/9906075 . 
C. Adler et al. [STAR Collaboration], Phys. Rev. C 66, 034904 (2 002) |arXiv:nucl-ex/ 0206001 . 

H. Satz, Quark Matter '02, Nucl. Phys. A 715, 3c (2003) arXiv:hep-ph/0209181 . 

H. Schlagheck [WA98 Collaboration], Quar k Matter '99, Nucl. Phys . A661 337c (1999). |arXiv:nucl-ex /9907005 . 
H. Sorge, Phys. Rev. Lett. 78, 2309 (1997) arXiv:nucl-th/9610026 . 

C. Pinkenburg et al. [E895 Collaboration], Phys. Rev. Lett. 83, 1295 (1999) [arXiv:nucl-ex/9903010l . 

D. Molnar and M. Gyulassy, Nucl. Phys. A 697, 495 (2002) [Erratum-ibid. A 703, 893 ( 2002) | |arXiv:nucl-th/0104073| . 
P. F. Kolb, P. Huovinen, U. W. Heinz, and H. Heiselberg, Phys. Lett. B 500, 232 (2001) arXiv:hep-ph/0012137 . 

M. Bleicher and H. Stocker, Phys. Lett. B 526, 309 (2002) arXiv:hep-ph/0006147 . 
Z. w. Lin and C M. Ko, Phys Rev. C 6 5, 034904 (2002) |arXiv:nuc l-th/0108039*r 

T. J. Humanic, |arXiv:nucl-th/0205053l . 

P. K. Sahu, N. Otuka, and A. Ohnishi, |arXiv:nuc l-th/0206 0Tol . 
D. Teaney, J. Lauret, and E. V. ShuryakT TarXiv:nucl-th/0110037| . 
J. D. Bjorken, Phys. Rev. D 27, 140 (1983). 

T. Hirano, Phys. Rev. Lett. 86, 2754 (2001) |arXiv:nucl-th/0004029| . 

H. Liu, S. Panitkin, and N. Xu, Phys. Rev. C 59, 348 (1999) arXiv:nucl-th/9807021 . 
S. Soff, S. A. Bass, M. Bleicher, H. Stocker, and W. Greiner, |arXiv:nucl-th/990306lf . 

L. V. Bravina, A. Faessler, C. Fuchs, and E E. Zabrodin, Phys. Rev. C 61, 064902 (2000) |arXiv:nu cl-th/0002045"l . 
M. Bleicher and J. Aichelin, arXiv:nucl-th/0205069 . 

N. S. Amelin, E. F. Staubo, L. P. Csernai, V. D. Toneev, K. K. Gudima, and D. Strottman, Phys. Rev. Lett. 67, 1523 
(1991). 

C. M. Hung and E. V. Shuryak, Phys. Rev. Lett. 75, 4003 (1995) arXiv:hep-ph/9412360 . 

D. H. Rischke, Y. P ursun, J. A. Maruhn, H. Stocker, and W. Greiner, Heavy Ion Phys. 1, 309 (1995) 
arXiv:nucl-th/9505014 . 

Y. B. Ivanov, E. G. Nikonov, W. Noerenberg, A. A. Shane nko, and V. D. Toneev, arXiv :nucl-th /00 1 1 004| . 

J. Brachmann et al., Phys. Rev. C 61, 024909 (2000) arXiv:nucl-th/ 9908010 ; J. Brachmann, A. Dumitru, H. Stocker, and 

W. Greiner, Eur. Phys. J. A 8, 549 (2000) |arXiv:nucl-th/9912014|. 

L. P. Csernai and D. Rohrich, Phys. Lett. B 458, 454 (1999) arXiv:nucl-th/9908034]. 

E. E. Zabrodin, C. Fuchs, L. V. Bravina, and A. Fassler, Phys. Rev. C 63, 034902 (2001) |arXiv:nuc l-th/0006056| . 
J. Gosset et al, Phys. Rev. Lett. 62, 1251 (1989). 



37 



[75] S. A. Bass, R. Mattiello, H. Stocker, W. Greiner, and C. Hartnack, Phys. Lett. B 302 (1993) 381. 
[76] B. A. Li, Nucl. Phys. A 570, 797 (1994). 

[77] J. Barrette et al. [E877 Collaboration], Phys. Rev. C 55, 1420 (1997) [arXiv:imcl-ex /9610006 . 



